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MATHEMATICAL STUDIES OF THE INTERACTION 
OF RESPIRATORY GASES WITH WHOLE BLOOD 
II. CO2 ABSORPTION. 


S. R. BERNARD* 


COMMITTEE ON MATHEMATICAL BIOLOGY 
UNIVERSITY OF CHICAGO** 


Some progress has been made on the problem of the interaction of re- 
spiratory gases with whole blood. A working mathematical model for the 
O2-CO2 interaction phenomena has been developed from mathematical 
studies of the data. The Edsall-Wyman (1958) model for CO2 absorption 
is improved upon in this paper by consolidating it with the O2 absorption 
model developed in paper I of this set (Bernard, S. R., Bull. Math. Bio- 
physics, 22, 391-415, 1960). This improved model assumed the effect of 
Oz on COg2 absorption is mediated through the electrical charge pos- 
sessed by the hemoglobin molecule, t.e., Og molecules bound to hemo- 
globin displace protons from the hemoglobin thereby increasing the nega- 
tive charge on the hemoglobin and at the same time increasing the acidity 
of the solution. The model is tested against the data. 


1. Introduction. In a previous report (Bernard, 1960), hereafter 
referred to as I, we developed and tested two models for oxygen 
saturation phenomena observed in studies with blood. From ex- 
tensive tests of one of these models against available experimental 
data we observed that we could fairly adequately represent the O2 
saturation curve over the pH range of ca. 6.5 < pH £8.88. We 
pointed out in I that the purpose in this research was to develop a 
working mathematical model for O2 and COz2 interaction with whole 
blood so that we could make some predictions with this model and 
from a comparison of predictions and measurements we determine 
the adequacy of the model. Now that we have a fairly adequate 
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model for Oz saturation phenomena we consider the problem of de- 
veloping a model for CO2 absorption phenomena. =. 

It is a well-known fact that the presence of oxyhemoglobin in 
blood suppresses the uptake of CO2. Here then is another test to 
which we can subject the model for O2 saturation, developed in I, 
but before we can consider this effect of oxyhemoglobin we need a 
model for COz absorption. 

Fortunately, a part of this problem may be considered as having 
been already worked out because Edsall and Wyman (1958) have de- 
veloped a working model for COz absorption. Their model consists 
of a simple solution of a salt, 0.02 molar sodium bicarbonate and a 
buffer, 0.025 molar methyl imidazole together with its hydrochloride 
salt also 0.025 molar. When the total CO2 contained in this solu- 
tion is calculated and plotted as a function of the partial pressure 
of CO2 gas to which the solution is exposed the resulting curve 
bears striking similarity to the curve for whole blood. In fact on 
plotting their theoretical curve together with the experimental meas- 
urements of Christianson e¢. a/. (1914) on whole blood the theoreti- 
cal curve superimposes on the measurements for the case where 
oxyhemoglobin is present in blood. (The buffer pk has to be 7.4 to 
make these superimpose.) Edsall and Wyman (1958) interpret the 
effect of oxyhemoglobin on the absorption curve in terms of the 
heme linked acid groups possessed by hemoglobin. These become 
more acid on oxygenation of the hemoglobin molecule. Their idea 
is that the pK of the buffer changes on oxygenation. It goes from 
8.0 or so, corresponding to the case where only reduced hemoglobin 
is present in the blood, to 7.4 corresponding to the case where only 
oxyhemoglobin is present. If one calculates and plots curves for 
these two different pX values the curves bear striking similarity to 
the curves observed in experiments with whole blood. 

This interpretation of the phenomenon is sufficient from a quali- 
tative standpoint but it is not satisfactory from a quantitative 
standpoint. It is not satisfactory because one has to arbitrarily 
change the pK value to get the theory to line up with the experi- 
mental findings. The Edsall-Wyman model can be rendered more 
adequate by postulating that the effect of oxyhemoglobin is on the 
electrical charge contained in the solution. We find that if we 
slightly modify their theory by adding to their model solution 0.002 
molar hemoglobin and assuming this to be charged, then on oxy- 
genation the hemoglobin becomes progressively more negative be- 
cause protons are displaced from the hemoglobin molecule on bind- 
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ing oxygen. The resulting theoretical CO2 absorption curves cal- 
culated for this case will coincide with the experimental curves. 

The number and the sign of the electrical charges possessed by 
the hemoglobin molecule, ¢.e., whether it is neutral or has four or 
only two positive charges, markedly affects the theoretical CO, ab- 
sorption curves. More COQz is contained by the model solution when 
uncombined hemoglobin is assumed to be positively charged and 
becomes neutral on oxygenation, while less can be contained when 
uncombined hemoglobin is assumed to be neutral and takes on a 
negative charge on oxygenation. 

The plan of approach in this paper is as follows. First, we dis- 
cuss the Fdsall-Wyman model. Then we present a very brief review 
of the findings we made in I and the models for O2 saturation phe- 
nomena. After that we modify the Edsall-Wyman model to incorpo- 
rate the O2 saturation model and we test the modified model against 
experimental data. Next we simplify the COe absorption equations 
for more convenient application to practical cases and for better 
interpretation. 

2. Mathematical Studies of CO2 Absorption. A. The EdsallWyman 
model for COzg absorption. Edsall and Wyman’s (1958) model solu- 
tion contains a salt, 0.020 M sodium bicarbonate, and a buffer, 
0.025 M methyl imidazole (Im) together with its hydrochloride 
(H*ImC1-) also 0.025 M. The chemical reactions for this system are 


BeOe Dae CO, eH. HCO; (1) 
HGOis 0 CO; (2) 
HIm* = H* + Im. (3) 


The terminology of Edsall and Wyman is employed here. Let p de- 
note partial pressure of CO2, T the total CO in the solution, 9’ 
the Henry’s law gas absorption coefficient for CO2, C the total 
imidazole concentration = 0.5 M, Kw, Ki, Ki and K the ionization 
constant for water, the first and second ionization constants for 
carbonic acid and the ionization constant for the buffer, respec- 
tively, and denote by / the hydrogen ion concentration. From the 
totality condition for the buffer 


C = [Im] + [HIm*] (4) 
and from the equation for electrical neutrality 
[HCO4] + 2[COs] + [OH7] + [CI] = [H*] + [Na*] + [HIm*] (5) 


and from the mass law expressions for the reactions in (1), (2) and 
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(3) we obtain the expression for p as a function of 4 which is 


p= CHAK + Abs h Ky/h + [No] eiClay (6) 
qd (Ki /h + 29Ki kg/h ) 


From the totality condition for CO2 


T = [COs] + [H2 COs] + [HCOs] + [C05] (7) 
and the mass law expressions for (1) and (2) we obtain 
T=p¢ (1+ Ki/h+ Kik3/h?). (8) 


Equation (6) is Edsall and Wyman’s (1958) equation (22) while (8) 
is their equation (10). 

Before consolidating this model with our model for Og saturation 
phenomena so that we may quantitatively predict the effect of O2 
on the CO2 absorption curve and test the same we give a brief re- 
view of the pertinent and relevant parts of I on the models for Oe 
saturation phenomena. 

B. Brief review of models for O2 saturation phenomena. In I we 
assumed that oxygen molecules and protons compete for binding 
sites on the hemoglobin molecule. From mathematical considera- 
tions of this model and in testing the model against the O2 satura- 
tion data we found that in effect, the following simple reactions 

ZX aja + 2X = ZXa; + (9) 
where 1<jS$2, and where Z and X denote molecules of hemo- 
globin respectively while H* denotes a proton, would adequately 
describe the interaction of O2 with solutions of hemoglobin as well 
as blood. This is called the Wyman-Allen (1951) hypothesis. We 
outline the derivation of the equation for saturation for we shall 
have use for it and the terminology later. 

Let lower case letters denote concentrations, ug; the concentra- 
tion of the ZX9;, Ka; the equilibrium constant for the jth reaction 

j 
in (9), Po; = 7m Kai the product of the equilibrium constants, 7, 
t= 


1 
the total concentration of hemoglobin 


Te =2+ Y) vay, 1 S782, (10) 
j 


T’, the total concentration of oxygen 


T, =@+2 D0 jusj, 1S 7S 2 (11) 


Z 
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and define 


S = (T, - x)/(4 Tz) (12) 


to be the fractional saturation (100S = percent saturation). From 
the mass law expressions for (9) 


Ua; = Poj(z/h®)?7, 1S j <2, (13) 
and from (10), (11) and (12) we obtain 
S= 1 > 93)P. 42; : 
eg ~. (29) Paj(a/h?)?i/ ¥* Poy (a/h4y2i (14) 
= j=0 


where Po =1. In testing this equation against available experi- 
mental data we found that 


P2/h = 0.8/2?, (15) 
P,/h? = 1/24 (16) 


where 250 is the oxygen tension for which S = 4. We tested this 


theory against measurements of saturation made on either solutions 
of hemoglobin or blood from four different species; sheep, ox, 
horse and human. In some cases the factor of 0.8 appearing in the 
numerator on the right hand side of (15) had to be changed to 0.96 
to obtain a better fit of the data. The agreement between (14) and 
the expressions in (15) and (16) with the experimental data was 
good enough for practical purposes. 

In I we also considered the case where oxygen molecules bind 
on the hemoglobin in four consecutive steps and only two protons 
are displaced from the molecule. We assume a random displace- 
ment mechanism for the proton displacement. This leads us to con- 
sider the case where more than one species of hemoglobin is pres- 
ent. At least six different species have to be postulated to ac- 
count for the release of two protons on binding four oxygen mole- 
cules. To illustrate this let the pair of numbers in parentheses 
correspond to the number of the oxygen molecule bound which dis- 
places a proton, ¢.g., (1,2) corresponds to the release of protons 
on binding the first and second O2 molecules. Hence we have 
proton displacement corresponding to the binding of the following 
pairs of O2 molecules, (1,2), (1,3), (1,4), (2,3), (2,4) and (38,4). 
This accounts for the six different species. One species releases 
protons on binding the (1,2) Oz molecules, another species of 
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hemoglobin releases protons on binding the (1,3) O2 molecules 
and so on. 

In the derivation of the equation for saturation for the case of 
more than one species, we assumed the reactions in which Og 
molecules were bound but did not release protons were equivalent 
and independent and characterized by a microscopic equilibrium 
constant denoted by a, while the reactions in which the O2 mole- 
cules displaced a proton were also equivalent and independent and 
characterized by a microscopic equilibrium constant denoted by 0. 
We made some simple and plausible assumptions for the weights to 
be assigned to these constants. The equation for saturation was 
found to be of the form 


% 
Sa1-9(1 + 2(B/h% + b2/B) ee (17) 


where § and y are functions of the a and 6 while z and A denote, 
as before, oxygen tension and hydrogen ion concentration, respec- 
tively. In application of (17) to O» saturation curves for human 
blood we found 6 = 107® and y=1075. The agreement between 
this theory and the measurements was excellent for O2 partial 
pressures generally greater than 20 mm of Hg. 

With this brief review of our findings in I we now return to the 
problem of modifying and improving the Edsall-Wyman (1958) model 
for CO absorption. In order to predict the effect of oxygen or 
oxyhemoglobin on the CO2 absorption curve in a quantitative man- 
ner, we postulate that this is mediated through the electrical 
charge possessed by the hemoglobin. Oxyhemoglobin becomes 
progressively more negative upon oxygenation because oxygen 
molecules displace protons from the hemoglobin molecule. Be- 
cause of this the equation for electrical neutrality, (5), must be 
modified. Now it makes a difference as to how we distribute the 
charge on the uncombined hemoglobin molecule, i.e., whether we 
assume uncombined hemoglobin has no charge or has one or two 
positive charges. This affects the shape of the CO, absorption 
curve as will be shown. 

C. An improved model for CO2 absorption. Assume uncombined 
hemoglobin is neutral. Denote by W the total negative charge con- 
tributed by oxyhemoglobin. W is given by 


W = HbOzZ + 2HbOz oe 


CO2 ABSORPTION i 
From (13) and (10) and (18) we obtain 


2 2 
W=T, (z arn) [5 Patent (19) 
j=0 J=0 


and from (14) we note that (19) can be written as 
W=27T,S8. (20) 


The equation for electrical neutrality for the case where oxyhemo- 
globin is present is 


[LHCO3] + 2[CO3] + [OH7] + [C17] + W = 
(H*] + (Na*] + [HIm*]. (21) 


From (5), (6), (20) and (21) and the mass law expressions for (1), 
(2) and (3) the expression for p as a function of A and 2, the O2 
partial pressure, can be derived and is found to be 


Ch/(K +h) +h-Kw/h + (Na*] - [Cl] -2T7,8 


p(x, h) = ? * Snecr 
Y (K{/h+2K{ K/h?) 


(22) 
Consider the case where reduced hemoglobin is present in the 
solution. Here, 2 in (22) is zero, hence $= 0 and then (22) is 
identical to (6). Thus p(0, 2) = p in (6). 
Consider the case where oxyhemoglobin is present in the solu- 
tion. Let 2, denote the O2 partial pressure in air required to satu- 
rate the hemoglobin. Then S = 1 and equation (22) becomes 


P(@1,h) =p—-2T72/(q’(Ki/h+2KiKg/h?)) (23) 
where again 7 is given by (6). 
In order to calculate T we have to employ equation (8) together 
with equation (22). We substitute p(a, A) for p in (8) to obtain 
T = p(z,h)q (1+ Ki/h+ Ki ké/h?). (24) 
Edsall and Wyman (1958) have made some calculations and we wish 
to make use of their data. In order to do so we use the following 
procedure. Specify a value for x and A and determine p(z, /) from 
(22). Then substitute this value for p(x, h) together with the corre- 
sponding value for / into (24) and determine 7. For the case where 
only reduced hemoglobin is present in the solution denote the total 
CO2 absorbed, for this case, by 7;. We proved above that p= 
p(0, 2) and hence we may write 
T, =p (1 + Ki/h+K{Kg/h?) (25) 
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where now 7 is given in (6). Edsall and Wyman have already cal- 
culated, tabulated, and plotted this function. Their values for 7, 
(they denote it by 7), p and A can be found in Table VI, p. 582 of 
their book. To calculate 7.x, the total CO2 contained in the solu- 
tion when oxyhemoglobin is present, we let x = az; in (24) and S=1 
in (22) and substitute for p(x1, 4) from (23) to obtain 


14K2/h+ Kiks/e (26) 
Ci/h + 2K KS/h? 


Tox = pg (1+ Kish + KiK3/h7) - 2T, ( 


From (25) and (26) we find 


- + ye 2 
dts Bae ath (27) 
K4/h + 2K4Ko/h 


Since 7, has already been calculated by Edsall and Wyman then 
from (27) T5x can be conveniently determined by inserting their 
values for 7, together with the corresponding value for h. 

A graph of T, and T,, plotted as a function of p and A appears 
in Figure 1, These curves were calculated with (8) and (27) and 
the Edsall-Wyman data. We assumed K =10~%* for the buffer, 


Voted pi ole 
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along with Edsail and Wyman and we let 1, = 0.002 M, the average 
concentration of hemoglobin in blood. Also plotted on this graph 
are the experimental measurements of Christianson, e¢ al. (1914) 
on the total CO2 contained in whole blood which was exposed to 
air, the case where oxyhemoglobin is present, and exposed to Hy 
gas, the case where only reduced hemoglobin is present. A large 
discrepancy between experimental measurements and theoretical 
predictions is noted. The reasons for this discrepancy are as fol- 
lows. Firstly, the theoretical curve for reduced hemoglobin, the 
Edsall-Wyman curve, fits rather closely the measurements for blood 
containing oxyhemoglobin. Secondly, the theoretical curve for the 
case where oxyhemoglobin is present was determined from the 
Edsall-Wyman curve. 

There are at least two ways to force better agreement with the 
experimental data. One way is to change the parameters in the 
Edsall-Wyman model such that their curve would superimpose on 
the measurements for the case where only reduced hemoglobin was 
present in blood. This maneuver would automatically bring the 
theoretical curve for oxyhemoglobin in closer agreement with the 
measurements. However, we do not know if this would make the 
curve coincide with the measurements. The difference between 
the theoretical curves can be seen from (27) to be due to the sec- 
ond term on the right hand side of (27). The factor 27, times the 
quotient of the polynomial functions of 1/A are responsible for the 
difference between the two curves. Remember the factor of 2 comes 
about from the assumption that saturated hemoglobin has two nega- 
tive charges. We point out had we assumed it had four negative 
charges, a situation which would arise if the four oxygen mole- 
cules displaced four protons rather than two protons, then the dif- 
ference between these theoretical curves would be much greater. 
Rather than undertake the formidable task of changing parameters 
in the Edsall-Wyman model—formidable because there are too 
many of them—we took an alternate tact. We determined what ef- 
fect positively charged hemoglobin would have on the model. Here 
we assume uncombined hemoglobin is a cationic acid. In deter- 
mining this effect we find out which parameter we have to change 
and we learn how to simplify the Edsall-Wyman equation for more 
practical application and more convenient interpretation of the 
data. 

When we assume uncombined hemoglobin has two positive charges 
while saturated hemoglobin is neutral the equation for W becomes 
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W =2 Hb + HbO, (28) 


= 22 + U2. 


Again, from the totality conditions for hemoglobin and oxygen and 
the mass law expressions we can write 


2 2 kas 
W=T, Y (2-f)Paj(a/h4)?1/ Yo Pa, (2/h*)”! (29) 
j=0 j=0 
where, again, Po = 1. The equation for electrical neutrality cor- 
responding to this case is 


[HCO3] + 21CO3] + [OH7) + [C17] = W+ [H*] + [Na*] + [HIm*] 


and when we substitute for W from (29) and solve for p(z,A) in a 
manner similar to that employed in obtaining (22) we find 


p(2,h) = ewe th) +haK,/h +Na-—Clo+ 
2 2 
T, )) @~ i) Paj(e/h*)7)/| | Pas (a/b)?! || + 
j=0 j=0 
[q’ (K4/h + 2K4K3/h?)). (30) 
When zx = 0 we have reduced hemoglobin and (30) can be written 
p(0,h) = p + 27, /(q’ (K/h + 2K4K9/h?)) (31) 


and when 2 =2;, the case where oxyhemoglobin is present, then 
(30) becomes 


Pp (@1,h) =p + Tz (P2/Pa) (h/e4)/(g’ (K/h + 
2K4K%/h*))= p,@1>>0. (32) 


In arriving at (32) we make use of the fact that only the terms in 
a, raised to the highest power in (31) are significant. Also we 
substitute for (P2/P4)(h/a7) from (15) and (16) to obtain 


(P2/P4)(h/e?) = 0.8 (e50/a1)* 


and since 259 for human blood is ca 25 while z; is ca 150, then 
this second term in (32) can be dropped. It contributes negligibly 
to the sum. Substituting into (8) from (32) for p we find 


Tox = pg’ (1+ K4/h + K4K3/h?) (33) 
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and from (31) and (8) we may write 


1.4 K%4/h + K4K3/h? 


ts = d Sr ste PY gs PS ae ES 
q’ (K4,/h + 2K4K3/h?) 


(34) 
where, again, 7, and 7.x stand for the total COz contained in the 
model solution when reduced and oxygenated hemoglobin is pres- 
ent, respectively. 

Figure 2 presents the graph of the theoretical curves plotted to- 
gether with the measurements for this case. As can be noted, 
theory and experiment are in close agreement. This tells us which 
of the parameters we have to change in the model if we assume un- 
combined hemoglobin is neutral and becomes negative on oxygena- 
tion, as we did previously, and observed a discrepancy in Figure 
1. All that is necessary is to add positive charge to the model 
solution in the amount of 2-(0.002) M and then the two curves in 
Figure 1 would be identical to those shown in Figure 2. Before we 
make use of this fact in simplifying the equations for more practi- 
cal application we consider another charge distribution. 
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Suppose 
+ = 
W = Hb" — HbO, (35) 
Stes Sie 
i.e., unsaturated hemoglobin has one positive charge, HbQ, is 
neutral and HbOs, completely saturated hemoglobin, has 1 nega- 
tive charge. For this case we may write 


2 2 

He Tye (Ls Pana heey (2 Pater (36) 
j=0 ime 

where Po = 1. When 2 =0, the case of reduced hemoglobin then 
W =+T7,. The curve 7, would then pass between the two theoreti- 
cal curves shown in Figure 2 passing below the measurements for 
blood containing reduced hemoglobin. But for the case of oxy- 
hemoglobin, let z = x; in (36) and then because fourth degree terms 
dominate, W becomes equal to-T,. The curve 7,,, corresponding 
to this case, would pass between the two theoretical curves in 
Figure 1, below the experimental measurements on blood contain- 
ing oxyhemoglobin. For some value of z in (36), W = 0 and then 
the curve 7 would be identical to the Edsall-Wyman curve in Fig- 
ures 1 and 2. 

For practical purposes the equations for p and T are too cun- 
bersome. We can make them a little more practical for application 
by making use of some of the findings we have made in the above 
plus some others which we have not mentioned. We could have 
presented this more practical model from the outset but we did not 
because we believe that it is more instructive to consider the 
separate cases in detail as we have done above. 

To obtain a set of simple equations for the CO, absorption phe 
nomena we assume that reduced hemoglobin is neutral and becomes 
negatively charged on oxygenation. This we learned will yield a 
curve which will undershoot the data (Figure 1) but this can be 
rectified easily by adding positive charge to the solution thereby 
raising the Edsall-Wyman curve, a ‘“‘pivot’’ curve, forcing it to 
agree with the measurements on blood containing reduced hemo- 
globin. Without too much loss in accuracy the term h — K/h, ap- 
pearing in the numerator of (22), together with the term q’2K/K¥/h?2, 
appearing in the denominator, can be omitted from (22) and the re- 
sulting equation for p(x, h) becomes 
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p(x,h) =(Ch/(K +h) + [Na*] -[Cl-] -2T, S]h/q’K?. (37) 


This expression is of the second degree in A, for the case where S 
is maintained at a constant value, and thus we can determine its 
root and express / as an explicit function of p and S. Let pq’K//K = 
p’, h/K =h’ and W’ = Na* —Cl7 -2T,S. Then from (37) we obtain 


p= PW 4+ Vip? - W244 (C +) 0} 
2(C +W’) 


(38) 


Now we can simplify the expression for 7, equation (24), by dropping 
the term K/K%/A?. It is negligible relative to the other terms in 
(24), for physiological pH values. Then we can substitute for A 
from (38) to obtain 


9 7? 7 7 
Ki (C +W’)/K 5) (39) 
p’ Ws Jip’ -W)2 + 4(C + 0} 

Figure 3 presents a graph of (39) for the case where S = 0, com- 


pletely reduced hemoglobin, and S = 1, completely saturated hemo- 
globin. We employed the following values in (39) to calculate T: 
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CG =0.05, K=10 "4, Na* =0.024, Cl” = 0.025 M, T, =0.002 M, 
gq’ = 1074-5 and K/ = 10761. These values are identical to those 
employed by Edsall and Wyman (1958) except for the Na* concen- 
tration. We increased this from 0.020 to 0.024 M to offset the addi- 
tion of negative charge contributed by oxygenated hemoglobin. As 
can be noted in Figure 3, the agreement is not as good as we ob- 
served in Figure 2 but it is believed to be close enough. 


The author is grateful to Professors N. Rashevsky and H. D. 
Landahl for their comments, suggestions, and criticisms. 
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ON THE MUTUAL INFLUENCE OF INDIVIDUALS IN A SOCIETY 


N. RASHEVSKY 


COMMITTEE ON MATHEMATICAL BIOLOGY 
THE UNIVERSITY OF CHICAGO 


On the basis of previously proposed mathematical models of social be- 
havior, the present note investigates the possibility of the control of be- 
havior remaining permanently in the hands of one class, if this class 
possesses sufficient means for influencing mass behavior. The conclu- 
sion is reached that, with the assumptions made, if the behavior imposed 
by the controlling class leads to sufficiently strong dissatisfaction, the 
control will pass to another class, no matter how strong the controlling 
power of the first. 


It is very frequently emphasized that in the so-called totalitarian 
countries not only is the behavior of every individual governed by 
the controlling group but the controlling group has in its possession 
such a powerful machinery of mass media influence as to make any 
counter-influence by dissenting individuals impossible. 

It has also been pointed out that such a situation, though in a 
different form, prevails in many so-called democracies. The be- 
havior of individuals is influenced through appropriate mass media 
not necessarily by the government but by other powerful groups 
which control the mass media. The ‘‘motivation research”’ gives to 
the controlling group a powerful tool to influence the masses with- 
out the latter even being aware of the influence (Packard, 1957). A 
change in the United States from a ‘‘society of opinions’’ of the 
past to a ‘‘society of masses’’ at present has been discussed by 
Wright Mills (1956). 

All this has found a reflection even in the fiction literature and 
led to such ‘‘utopias’’ as pictured in George Orwell’s 1984 or 
Aldous Huxley’s Brave New World. The situation depicted in those 
novels implies actually a complete standstill of any social change 
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due to a continued control by one class which holds all the con- 
trolling apparatus in its hands. 

In previous publications (Rashevsky, 1947, 1959) we developed a 
mathematical theory of the influence exerted by a smaller class of 
‘‘actives’’ on a larger class of ‘‘passives.’’ We found that various 
factors may cause a shift which passes control over the behavior 
of the majority from one group to another. In all such cases, how- 
ever, we considered that the ‘‘coefficients of influence’’ of the two 
groups which compete for control are of the same order of magnitude. 

The purpose of this note is to investigate whether, on the basis 
of the mathematical models developed by us previously, the outlook 
mentioned in the introductory paragraphs above is justified. 

First we shall consider the simple linear interaction which is a 
first approximation to a more elaborate model (Rashevsky, 1960). 

The condition for the control of behavior by class I of Xo active 
individuals which competes with class II of Yo actives is 


Xo: SENG eee (1) 
ao + @1 ao + Qi 
where do is the coefficient of influence per individual of the first 
class; Co, the coefficient of influence of the second class; and ay, 
the corresponding quantity for the passives. 
The inequality 


) ne ES AP, 
, Co + @ ee tas : (2) 
insures control by class II. 
Changes in control occur as a result of either changes of Xo or 
Yo, or of ao and co, or of all those quantities. Suppose we have 
originally inequality (1). Suppose further that 


@o — @1 > No(co + @1). (3) 


Then even if Xo = 1, that is if the first class consists of only 
one individual, inequality (2) can never become satisfied because 
it would imply Yo >No. If, as we always assume, a1 << Go ; 
@1 << 0, condition (8) requires that @o be at least No times co. 
With No of the order of magnitude of 107-108, it may be questioned 
whether such a tremendous influence can be held by even the most 
strongly organized group. In principle, however, the situation is 
possible, especially considering the effects of modern mass com- 
munication media and the practical impossibility for individuals 
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who are not members of the controlling group to operate those media 
on a scale even remotely comparable to that used by the controlling 
group. 

The situation remains very much the same in the more general 
case (Rashevsky, 1960) so far as the direct effect of the active in- 
dividuals is concerned. An important difference, however, occurs 
because of the nonuniformity of passives. 

As we have seen (Rashevsky, 1960, p. 70), the conditions for 
behavior F , (control by the first class) are in this case, using the 
same notations as before, 


abo + CoY 
bg el aR 
ao 


2 VseNs —a)3 
BO sts 3ab3BNo (5) 


The condition for control by the second class is 


(4) 


where 


Yoo abo + AnXo (6) 
Co 

ence, if (4) holds and if @ > Noco, inequality (6) will never 
occur. 

ilowever, in the present case expression (5) holds only when 
N (¢@) the distribution function of the individual tendencies of the 
passives toward behavior Ff; or Ro, is symmetric with respect to 
¢ =0 and has a maximum at that point. As we have seen else- 
where (Rashevsky, 1959, chap. xxii), N (¢) may change as a result 
of ‘‘social learning.’’ If the shift is sufficiently large in the di- 
rection of negative values of ¢, then bo becomes negative and de- 
creases to -~ as 6 —»-o. Hence, no matter how large co, a 
sufficiently negative bo can always make (4) satisfied. 

As we have seen (Rashevsky, 1959, chap. xxii), the shift toward 
negative values of ¢ occurs when the activity 7, leads to dissatis- 
faction of the individuals. This shift, however, does not proceed 
indefinitely but tends asymptotically toward a finite value, at least 
under the assumptions which we made in Joc. cit. This asymptotic 
value is determined by the different parameters which enter into the 
theory. Therefore bo also tends to a negative asymptotic value. 
‘If this value is such that |a@bo|>aoXo, then inequality (6) will 
eventually become satisfied, and the control will pass to the other 


class. 
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Thus the non-uniformity of the passives with respect to ¢ and the 
shift of ¢ as a result of dissatisfaction with the activity imposed 
by the controlling class will eventually lead to a change of control, 
no matter how strong the latter. If, however, the behavior imposed 
by the controlling class is satisfactory and pleasant, then no 


change will occur. 


This work was aided by a grant from the Dr. Wallace C. and 
Clara A. Abbott Memorial Fund of the University of Chicago. 
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Traffic in one direction on a multilane highway is considered, and a 
general expression for the number of cars which pass a car travelling at 
a given velocity, as well as the number of cars which the given car 
passes, is derived for the case when the speeds of different cars are 
distributed in some arbitrary manner. Closed expressions are derived 
and discussed for a rectangular distribution. Each passing by another 
car or of another car is considered as a distracting stimulus which af- 
fects the reaction times of the driver. Using previously derived expres- 
Sions for the safe speed as a function of reaction times, expressions for 
the safe average speed are derived, in terms of the volume of traffic and 
of the spread of the distribution of speeds. 


In previous papers (Rashevsky, 1959a, 1959b, 1960; hereafter 
referred to as I, II and III) we discussed some biomathematical as- 
pects of driving on an empty highway. We found that conditions of 
safe driving require definite relations between some purely me- 
chanical and geometric parameters (such as width of the highway 
lane, roughness of pavement, etc.) and purely biological ones (such 
as reaction times). In the second of the two papers we discussed 
in particular the effect upon reaction time of various ‘‘distracting’’ 
stimuli, which may be either inherent in the driving conditions, 
such as unevenness of the pavement, or attributes of the highway 
surroundings, such as signs on the road, etc. 

In the present paper we shall consider driving in heavy traffic 
on a multilane superhighway, in which every car is frequently both 
passed by, and itself passes, other cars. As a rule traffic is not 
too heavy on multilane superhighways. However, situations dis- 
cussed here do occasionally occur. They seem to be rather the 
rule on the expressways in the neighborhood of Detroit or on the 
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Queen Elizabeth Highway between Toronto and Hamilton, Canada, 
as well as on some California freeways. This situation also cre- 
ates distracting stimuli, and we shall apply to it the method used 
previously (Rashevsky, 1959b). The readers’ familiarity with both 


preceding papers is assumed. 
Consider first for simplicity two lanes in one direction. Let the 


total traffic flow of the road in one direction be denoted by C cars 
per unit time. Let the average speed be 0. For a random distribu- 
tion of cars along the highway this requires that the average num- 
ber N of cars per sufficiently large stretch of highway be 


ye (1) 


or 


U= 


(2) 


ZI 


For each lane the number of cars per unit length is = N. Hence 
neglecting the finite size of the car, the average distance between 
cars is 

Fives 
l=, 
= (3) 


If the average length Jo of the car is measured in the same units 
as the length of the highway, we have always 4Nlo < 1. The total 
free length of the lane not occupied by cars is 1 - 4Nlo; hence the 
average distance between the cars is 


— = lo. (4) 


Considerations of safe stopping require that the speed v satisfy 
equation (15) of I, which in the present situation is 


G0 * 407s a hae (5) 


Here a is the constant which determines the stopping distance p of 
the car in terms of » by the equation 


p= av?, (6) 


and 7; is the reaction time to the stopping of the car in front. The 
positive root of equation (5) is 


wee —T1 + Vier nn 4al (7) 


Qa ; 
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Introducing equation (4) into equation (7) we find: 


ty + Vtr ~- St —4alo 


v= ae (8) 


This equation expresses the maximal safe average speed in terms 
of the density N of the cars, the size of the car lg, the mechanical 
constant a, and the biological constant T;. But, according to 
equation (1), the density N is itself a function of v. Introducing 
(4) into equation (5), we find 


= 9\ 
w+ (4 2)o+t=0. (9) 
If 
2 
Ce — 
1” (10) 
so that 
se $ +0, (11) 


then equation (9) has no positive real roots. For sufficiently small 
values of Jo, both roots are real and negative. For large values of 
lo, they are complex. Hence we must have 


2 
C< => 
= (12) 
But when (12) is satisfied, then (9) has two positive roots, a re- 
sult which at first sounds paradoxical. The paradox is resolved 
by remembering that the condition of safe stopping is given not by 

equation (5) but by 
ce ee 1g 0 (13) 


[cf. expression (15) of I]. Therefore instead of (9) we should write 
av” Fe (- 2a lo <0 (14) 


as conditions of safe stopping. But, if v; and v, > vs denote the 
two positive roots for the equality sign of (14) when (12) is satis - 
fied, then the inequality in (14) is satisfied only for 


Deane < V1. (15) 


In other words, a car going at a speed 7 < v, will not satisfy 
(14) and therefore will not satisfy (13). It will not be able to stop 
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in time because for such low values of v for a given total traffic 
flow C, the distance between the cars will become too small even 
for velocities 0 < v;. Hence there is both a lower and an upper 
safe speed for a given total traffic flow C, as long as that total 
traffic flow satisfies (12). If it does not satisfy (12) no safe speed 
exists. 

The conclusion about the existence of a lower limit for » may 
seem absurd because we could argue that at » = 0, everything 
would be safe. But for fixed C, the condition v = 0 would require, 
according to (1), N = « which, according to (4), would require a 
negative distance 7 between the cars, an obvious physical im- 
possibility. 

Such dense traffic, where the lower limit v; is approached (from 
above), is not likely to occur; therefore we may safely disregard 
the smaller root of (9) and take for the safe average speed the 
positive root 


(16) 


For Jo = 0, the smaller root of (9) becomes zero, while the larger 
one, given by (16), becomes 
— 2 rie 
oie Or aac (17) 
Now, however, we must remember that, as shown in II, the reac- 
tion time T; is also an increasing function of the velocity. On an 
empty highway this increase is due, as we have seen, to uneven- 
ness of the pavement, to wind noise, and to distracting stimuli of 
other types. If on a two-lane highway the speed is distributed 
around the average, then to the distracting stimuli studied in II are 
added the stimuli caused by the necessity of passing another car 
as well as by being passed by another car. 
Let the traffic speed be distributed around ¥ and let the fraction 
of cars with a speed between v and v + dv be given by 


f(v) dv (18) 
with 


i f(v) dv =1. (19) 
0 


saat I lad 
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If the speed of all cars were the same, in other words if the 
standard deviation of f(v) were zero, then there would be no pass- 
ings. Let there be n cars per unit length all going at a speed v, 
and let one car go at a speed vw > v. That car would overtake all 
the other cars. How many cars will it overtake per unit time? That 
number will be the same as if all cars were stationary and the one 
we consider were moving at a speed wv’ —v(>0). If there are n 
cars per unit length of highway, then the number of cars overtaken 
by the ‘‘exceptional’’ car would be per unit time 


n(v’ — v). (20) 


Hence, if the cars are distributed along the highway randomly as 
regards speed, a car going at a speed v’ will overtake per unit 
time a number of cars going at a speed between v and v + du(v < v’) 
given by: 

(v’ - ») Nf(o) de, (21) 
where, as before, N is the total number of cars per unit length of 
highway. 

Hence the total number of cars passed by a car moving with a 
speed v is given by: 


Py (v’)=N fe — v)f(v) dv. (22) 
0 


By a similar argument we see that the total number of cars which 
will pass a car going at a speed wv’ is 


Pa(v’)=N is (v — v’) f(v) dv. (23) 


As an illustration only, we shall derive expressions for P; (v’) 
and P(v’) for the simple unrealistic case in which f(v) is a rec- 
tangular function. In other words, we assume that for v1 < v< v2, 
we have f(v) = c(c = constant), while for v < v; and v > va, f(v) = 0. 
Requirement (19) then gives 


pa ae (24) 
Vg —- V1 


so that for v1 < v < va, we have 


f(r) = ’ (25) 


V1, — v2 
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The average speed v is then given by 


> = + (26) 


Introducing (25) into (22) and (23) and remembering that in (22) 
the lower limit of integration is now v1, while in (23) the upper 
limit is v2, we find 


P,(v’) = N foe (27) 
and 
Pa (v’) = N a= vy (28) 
2 (v2 - v1) 
We always have 
V1 <U' <2. (29) 


Therefore in (27) v’ — v1 < v2 — v1. When the spread v2 — v1 tends 
to zero, P;(v’) also tends to zero. The same thing holds for P2(v). 
For zero spread, when v = v1 = v2, there are no passings. 

When v’ = 01, Pi(v’) is zero. The car driven at that velocity 
does not pass any other cars. But for that value of wv’ the value of 
P2(v’) is largest. It is equal to 


P2 (v1) = Aiea = 2), (30) 


Similarly, for v’ = v2, Pi (v’) is largest, having the value 


Py (va) = M@a= 0), (31) 


while P; (v2) is zero. 
The sum P(v) = P; (v’) — P2(v’) gives the total number of cars 
either passed by or which pass the given car. It is equal to 


P (v') = won)? + 0" = 91)? 
2 (2 — v1) 


This quantity has a minimum for 


(32) 


v = ere (33) 


Hence a car going at the average speed has a minimum of the total 
number of distracting stimuli due to either passing or being passed. 
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It is convenient to express in (27) and (28) v; and vg in terms of 
v and of the spread 


QA = %. — v1. (34) 
From equations (26) and (34) we have 
V1, =u-A; Wm =04 A, (35) 


Now equations (27) and (28) become 


P, (0) = N(w = 7+ A)? 


r N(w +A —v’)? 
Py (v7) = N@*E =e) (37) 
while P(v’) [eq. (32)] becomes 
—2 “Eealteee 72 
P(v’)=N ia St alas (38) 


2A 


Expression (38) shows immediately that P(v’) has a minimum for 
U =. 

The reaction time T; to the stopping of the car in front is given 
by an expression similar to expression (2) of II for the reaction 
time to the deviation of the car from the straight course. The 
stimuli and threshold being now different, we mark them by the 


same symbols as in II but with a subscript ‘‘1.’’ Thus 


ig ee (39) 


This reaction time T; is affected by the distracting stimuli dis- 
cussed in II in the same manner as T. The role of T; is, however, 
somewhat less critical. If T; is too large, a driver will not be able 
to stop, at a given speed, his car in case of a sudden stopping of 
the car ahead of him, and a collision will occur. But zf the car in 
front does not happen to stop suddenly, then the driving may pro- 
ceed safely for an indefinite time, even though 7; is too large. On 
the other hand, if tT becomes so large that the right side of equa- 
tion (1) of II becomes smaller than the left side, the driver will 
lose control of his car regardless of what other cars do and even 
regardless of whether there are any other cars on the road. 

Both the passing of another car and the passing by another car 
are distracting stimuli of the nature discussed in II. They are, 
however, in general not equal. Passing another car is likely to be 
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a much greater ‘‘distraction’’ from relaxed driving than being 
passed by another car, when the distraction may be due merely to 
the noise of the other car. Therefore, if we write for the contribu- 
tions of both stimuli to R in equation (2) of II for t the effects of 
passing and of being passed €P; (v’) and CP2 (v’) correspondingly, 
we shall have €>¢. The other contributions discussed in II re- 


main and we have for T: 


ti Ey 
= i log r= Ro Sd, yu’? ie. nv’ - EP, (v’) Ey CP» (v’) . (40) 

In III we have generalized the discussion of II for the case in 
which the angle @ is not constant. This more realistic case leads 
not to a sharp critical speed, beyond which safe driving is impos- 
sible, but rather to an expression for a probability of an accident 
a function of speed. The expression reduces to that obtained in I 
and II, when the standard deviation of 9 is very small. For simpli- 
city we shall in the following consider this limiting case. The re- 
sults can be readily generalized for the more general and realistic 
case. 

Introducing (40) into equation (1) of IJ, and putting w’ = v, we ob- 
tain, as before, an expression for the maximum safe average speed 
on a multilane one-way highway as a function of the average speed 
and of the spread of speed. The equation thus obtained is in gen- 
eral of fourth degree and its closed solution is rather cumbersome. 
As an illustration we shall consider the case where the term yv’5, 
which is due to the unevenness of the road,can be neglected. This 
is likely to be the case on a new well-maintained superhighway. 
We shall also make the restrictive assumption used in II, that 
E, >> Ry, so that we can expand the logarithm and preserve only 
the linear terms as we did to obtain equation (16) of II. Then, 
from equations (36), (87), and (40), putting 


f+ C=, (41) 


we find 


_ . So + 4nv + NxA 
4kE (42) 


Introducing this into equation (1) of II, remembering (2) of this 
paper and using (18) of II, we find after rearrangements: 


4nd? + 4R ov + (XCA — 4kAE) = 0. (43) 


It must be kept in mind that the parameters &, Ro, n and E vary 
from driver to driver even if they all drive at the same speed. The 
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constant A varies from car to car. Hence in determining the aver- 
age v by means of equation (42), we must also take average values 
for the above parameters. Therefore in all the following equations 
we shall mark them with bars. 


If 
xCA - 4kAE < 0, (44) 
or 
4kAEF 
ie 
eAT Ye (45) 


then equation (43) has one positive (and one negative) root. If in- 
equality (45) is not satisfied, then equation (48) has only negative 
roots, which means that in that case the only safe speed is zero. 
If A =0, then the safe speed, given by (41), becomes independ- 
ent of C, and C could be made as large as we wish. This is due 
to the fact that we considered the effect of passing only on T, the 
reaction time for deviation from the straight course. We must, how- 
ever, also consider the effect upon 1, the reaction time to the 
stopping of the car just ahead. The maximal safe speed with re- 
gard to sufficiently rapid stopping is now obtained by introducing 
the appropriate expression for T; into (16). The reaction time 71 
is represented as a function of the other parameters by an expres- 
sion similar to the one which represents T, in other words by an 
expression obtained from (40) by merely adding the subscript ‘‘1”’ 
to all the parameters. Using such an expression, omitting again 
the cubic term yv3, and introducing it into (16), we find a fourth 
degree equation for the determinants of v. A simpler expression 
is obtained by using the approximation (17). We then have 
8k,;E,;-4CRo, — 


4(ak,E, + 71)0? - v+x1CA=0. (46) 


C 
If 
8k, E- 40R; 
cet ef)’ (47) 
C 
that is, if 
G > 2he (48) 
Ro 


then equation (46) has no positive roots. If the inverse inequality 
holds, then there are two positive roots, and as before the maximal 
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safe speed is given by the larger one. The smaller one again 
represents the ‘‘minimal’’ safe speed (see p. 21). Putting 


4(ak,E, + m1) =D, (49) 


we have for the maximal speed, from (45): 


ki Ey - 20R Fan oe, Ne 
ee (a8 20Ro, (4k Bi ) aed ah 
D C C 


In order for a safe speed to be possible, v must not only be posi- 
tive but alsoreal. Therefore, in addition to the condition (48), we 
also have 

4%, E, -2CRo,\? . 

(ee Si) See ae 


or 
1 eee: ~ sie 


Puc C (52) 


For constant C, v increases with decreasing A, as seen directly 
from (50). For A = 0, the total traffic flow C still must not exceed 
2k1Ei/Ro,, according to (48). Thus, as we anticipated, in this 
case even for A = 0, the safe speed v depends on the total traffic 
flow C. As C increases, v in this case decreases, becoming zero 
for C —» 0. When A #0, the safe speed also decreases with in- 
creasing C. 

It will in general depend on the numerical values of 7, 71, Ro, 
Ro,, E, E1 etc. which of the two safe speeds, that given by (48) or 
that given by (50), is larger. The smaller of the two values is to 
be taken as the actual highest admissible average speed. 

Equations (438) and (50) can be verified experimentally by use of 
devices which simulate driving conditions. The experiments will 
not be simple to set up, but they are more than likely to pay large 
dividends by contributing to our knowledge of the conditions of 
safe driving. 

Both according to (43) and to (50), the safest thing to do is to 
make A = 0, that is to maintain the same speed for all cars. It 
may be suggested that with the further increase of highly traveled 
highways, an installation of governors on engines, which would 
maintain a prescribed speed, would be made compulsory. Or else, 
from equations (43) and (50), we may calculate the admissible A 


for desired values of 7 and C and attempt to maintain this value of 
A by some means. 


i) Tih ee 
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Instead of putting v’ = v in (38) and in the corresponding expres- 
sion for T1;, we may, as mentioned above (p. 25), find a relation 
between v’ and v and A which will indicate what percentage of 
cars, for a given spread A and a given average speed 2, are driv- 
ing beyond safe speeds and thus are most likely to meet with ac- 
cidents. Such calculations could be compared with available data. 

The above conclusions hold only for the assumed ‘‘rectangular’’ 
distribution of speeds. Except for more complicated calculations, 
a similar study can be made for any distribution function. The ex- 
pressions will become more complicated, but the general character 
of the relation between v and the standard deviation A is not likely 
to be much affected. This is one of many problems for future 
studies. 

As in II, we can also introduce not only the instantaneous ef- 
fects, but integral fatigue effects. This will be made the object 
of a separate study. 

In conclusion we may remark that similar considerations as above 
can be applied to two-lane two-way highways. If the total traffic 
flows, which may be different in each direction at a given time, 
are given, then the number of cars which a given car meets per unit 
time is readily given. This number may then be considered as a 
distracting stimulus, affecting tT. The expressions thus obtained 
will give a dependence between the safe maximal average speed 
and the total flows in each direction. 


This work was aided by a grant from the Dr. Wallace C. and 
Clara A. Abbott Memorial Fund of the University of Chicago. 
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This report describes a method by which the concepts of the lead vec- 
tor and the lead field may be extended to include equivalent cardiac sin- 
gularities above the first order. The theoretical treatment results in 
generalization of the well-known lead vector principle into a lead tensor 
_ apparatus. The specific electrocardiographic contribution of an mth order 

multipole is found to be the product of an mth rank contravariant (‘‘heart’’) 
tensor by a covariant (‘‘lead’’) tensor of the same rank, and with identi- 
cal indices. Methods are also described whereby lead tensor components 
may be cast into a form which relates them directly to the more usual, 
spherical harmonic notation of multipolar current generators. Some atten- 
tion is devoted to the possibility of applying lead tensors to biomedical 
problems, particularly to the prospect of quantitatively evaluating the 
significant multipolar components of the human heart. 


It is common practice in analyzing the relationships between the 
electromotive forces of the heart and various electrocardiographic 
lead connections to represent the electrical behavior of the heart 
as one or more current dipoles. The artifice of the equivalent car- 
diac dipole is at best a first approximation, especially when it is 
limited to the case of a single, fixed-location dipole. More re- 
cently, idealization of the heart as a fixed-location current singu- 
larity has been extended by the addition of higher order terms to 
the dipolar term (Yeh and Martinek, 1957, 1959; Geselowitz, 1960). 
This generalized form of the equivalent singularity concept ap- 
pears to be well worth investigating since it can probably be made 
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to account in a strict quantitative manner for the deflections of any 
conventional electrocardiographic lead. 

In attempting to determine true heart-lead relationships, the bio- 
medical scientist is confronted with a two-fold problem. In es- 
sence, he wishes ultimately a) to determine the intrinsic electrical 
properties of the heart itself, and b) to evaluate the relationship of 
the heart to the irregularly shaped volume conductor within which it 
is normally contained. The problem appears hopelessly complex 
unless at least one of its two aspects can be treated independently 
of the other. It is essentially for this reason that the heart is 
idealized into some sort of equivalent electrical generator, thus 
permitting the volume conductor properties of the body to be de- 
termined relative to an arbitrarily selected type of heart model. 

The intrinsic properties of electrocardiographic connections are 
well known in relation to the equivalent cardiac dipole, and are 
best expressed as lead vectors (Burger and van Milaan, 1946-48). 
Within the framework of the equivalent cardiac dipole approxima- 
tion, an electrocardiographic deflection is defined as the scalar 
product of a heart vector and a lead vector. Although the lead vec- 
tor concept has added greatly to the accuracy of electrocardiographic 
analysis, it is of itself only a first order approximation since it 
fails to account for the electrocardiographic contributions of equiva- 
lent singularity components above the first order. The present com- 
munication describes our attempts to obtain more accurate charac- 
terization of lead properties by generalizing the lead vector concept. 
As might be anticipated, the electrocardiographic contribution of 
each higher order component proves to be the scalar product of an 
equivalent heart tensor by a lead tensor. 

The desired theoretical development follows rather directly from 
lead field considerations. By definition, the McFee-Johnston lead 
field is the electrical field which is produced in the body when one 
unit of direct electrical current is causéd to flow from the positive 
to the negative terminal of an electrocardiographic connection 
(Johnston, 1957). Because of the Helmholtz principle of reciprocity 
the gradient of the lead field is identical to the lead vector (Brody 
and Romans, 1953). Accordingly, the electrocardiographic contribu- 
tion of a dipole of moment, M, is M-F, where F is the lead field 
gradient. 

In rectangular cartesian coordinates, X!, X2, X 8, with orthonormal 
base vectors, a', a?, a3 (which are identical to @1, G2, as), the 
lead vector, F, and its successive differentials are of the form, 
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Oa? 
2 
dF =a! aN - dal 
Ox'da! 
(1a) 
dx? Oxi dak 
, 4 : 
BiB at ay dei dx* dz!, etc. 


? 


dxidax! dx* dx! 
dai = i ds® 

dx* = )k ds (8) (1b) 
dz! = )!'ds, ete. 


where wu is the lead field potential, and the A symbols are the di- 
rection cosines of the differentiation pathways. 

By definition, an nth order multipole is formed by bringing two 
mn —1 order multipoles of equal magnitude and opposite sign into 
infinite proximity. The moment of the multipole thus produced is 
defined as H yn) = nHcn_1) ds” where ds“ is the interpolar sepa- 
ration of the two lower order multipoles. In order for H:,) to have a 
finite value, H:,-1;, must become infinitely large as the two n - 1 
singularities are caused to approach each other (Stratton, 1941). 

Beginning with first order singularities, consider a dipole of mo- 
ment, M = a;M! =a4;H,1) “PA'. The electrocardiographic contribu- 
tion of the dipole is Va) =M-F = MF; = Hay PA‘ du/de'. Pro- 
ceeding to higher order singularities, Vi2) = — M«dF, V3) = M-d?F, 


Vay = -M-d'F,....... 220) Vay =(-1)"7'M-d"-'F, from which 
we obtain by the substitution of equations 1, 
1, 07u 
V =a Cyt (2)yI é 4 ds?) 
4 se dx' dx! 
_ 4 yi @s 074 . 
dx ' da) 
A) (yi (2)yi 0° u (2) 
= — Ai ye —_—_ 
aay mia dx* dai dak 


1 oe H (a) (1)yi @)yj (3) k (ADI 
Vay se 4! 


dtu é 


34 D. A. BRODY, J. C. BRADSHAW, AND J. W. EVANS 


In the preceding equations, and in those to follow, the summation 


convention is observed for 2, 7,4, ..+++-+ee> = 1,2, 3. 
Two sets of symbols may be now defined as follows: 
F; = du/dx' 


Bige = Otu/ aa Aa 


(3a) 
Fi jx = xy Ou/da' dai dk 
Fijkt = fa d4u/ dx dxi dx* dx', etc., 
Mt = Hay AN 
ee (3b) 
Miik — H i) Atik 
miikl — H,, AUF! etc., where 
Ai = (yi 
tj, BM) A 
. : es (3c) 
Atik oa ae (3)\k Au 
Aiikl | _~ A)! At ete, 


The F symbols are recognized as covariant tensors, and the M 
symbols as contravariant tensors. It is evident that the electro- 
cardiographic contribution of an nth order singularity is the product 
of two tensors of the nth rank. For instance V4) = M‘*! x Fy 57. 
Since the indices of the parent tensors are identical, their product 
is a scalar. In order to generalize the nomenclature of the equiva- 
lent dipole concept, we refer to the M symbols as heart tensors, 
and the F symbols as lead tensors. Indeed, the first rank tensors, 
M‘ and F;, are recognized as heart and lead vectors, respectively. 
The A symbols, which are formed by successive outer multiplica- 
tion of direction cosines, represent unit strength multipoles in 
tensor form. 

Because the order of partial differentiation is immaterial, there 
are +(n + 2)(n + 1), rather than 3", ap Sreany distinct components 
of the nth rank lead tensor. Of these, ¢n(n — 1) are redundant be- 
cause the lead field is a conservative nbtsatial field, leaving 2n +1 
components sufficient for lead tensor quantification. 
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Although the tensor apparatus provides an extremely compact and 
powerful means of expressing the heart-lead relationships of higher 
order singularities, a more explicit type of treatment is required to 
cast it into specifically usable form. One possible way of meeting 
this requirement would be to express all the spherical harmonic 
components of an nth order singularity in the form of corresponding 
M symbol components. However, we have found it more rewarding 
to gather the lead tensor components into subsidiary groups which 
bear a clear and direct relationship to the spherical harmonic nota- 
tion of potential functions. 

The desired regrouping is determined by analyzing lead fields of 
the form, 


cos md 


1 Wale ce i ee 


1+5,° (n+m)! sin m¢ Pa® (u) (4) 


where ru =r cos 0=z2=2!', r sin 0 cos 6=y=~2?7, r sin @ sin 
@= 2-2", and the P(«) are associated Legendre functions. Be- 
cause equations 4 are of the nth degree in 2, y, and 2 (cf. Byerly, 
1893, pp. 206-7), the associated lead tensors have the following 
special properties: (i) components above the nth rank vanish, (ii) 
components of the nth rank are constants (or zero), (iii) existing 
components of less than the nth rank are position-variable, but 
vanish at the coordinate system origin. Because components of the 
nth rank are invariant, equation 4 lead fields of the nth degree are 
‘‘ideal’’ with respect to singularities of like order in much the same 
sense that uniform (i.e., first degree) lead fields are ideal with re- 
spect to the dipole moment of source-sink current distributions 
(Brody, 1957). 

Tables I, II, and III list, respectively, the F; of the ideal first 
degree lead field, the F;; of the second degree field, and the Fijx 
of the third degree field. From these, and corresponding higher 
rank components, we form the scalar products, MF, Heath a 
Ati‘ F-., etc. In the second rank case, for example, the following 
products are obtained: 


m cosine series sine series 

0 12 — 7 (81 B2 + yr y2) 

1 4 (41 B2 + 21) x (1 y2 + A2y1) 
1 4 (81 B2 - y1ya) 4 (Biya + Bayt) 
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TABLE I 


FIRST RANK TENSOR COMPONENTS OF 
FIRST DEGREE IDEAL LEAD FIELDS 


(cf. equation 4) 


series 


TABLE I 


SECOND RANK TENSOR COMPONENTS OF 
SECOND DEGREE IDEAL LEAD FIELDS 


(cf. equation 4) 


where G1, B1, yi = + “PA*, and Gg, Bo, yo = + “A. These prod- 
ucts are recognized as the coefficients, a2, and bem, of the in- 
dividual zonal and tesseral harmonic components of a unit strength 
quadripole (Stratton, 1941). A similar relationship exists between 
all other equation 4 lead fields and the harmonic component co- 
efficients, @,, and b;m, of corresponding order and degree.* 
Because of the specific characteristics of equation 4 lead fields, 
as illustrated above, Tables I, II, III, etc. serve as guides for the 
reduction of lead tensor data derived from other lead fields, hypo- 
thetical or experimental. In essence, all that is required is to sort 


*The implication is not intended that the notation, a@,, and b 
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and tabulate tensor components of a given rank into rows and col- 
umns in such a manner that division of each row by a weighting 
factor, Pam OF Ynm, Will reduce the tabulated data to the ideal table 
of corresponding rank. In principle, then, an electrocardiographic 
deflection at any given instant can be expressed as 


n=co m=n 


Y= : Pn an + ye YE Go Qnm> Ynm bam) (5) 


n=1 n=1 m=1 


When limited to terms of the first degree, equation 5 reduces to 
V = pia, + 11011 + Gi1011, which is essentially identical to the 
linear transformation form in which Burger (1946) expressed the 
case of the equivalent cardiac dipole. Conversely, when terms of 
all order and degree are included, equation 5 becomes a generali- 
zation of Burger’s equation. Assuming that the heart can be con- 
tained in a spherical surface which does not intersect the body 
surface, equation 5 can probably be derived from the principle of 
superposition. However, we have elected to present a derivation 
based on the principle of reciprocity because of the prospects for 
experimental application which are implicit in this approach. 

Because the lead fields given by equation 4 are specifically re- 
lated to individual harmonic components, they have certain theo- 
retical implications which are worth mentioning briefly. It will be 
recalled that compound electrocardiographic connections which pro- 
duce uniform lead fields (i.e., lead fields which are ‘‘ideal’’ in the 
first degree) respond only to dipolar singularities. Accordingly, 
this property permits the dipole moment of source-sink distributions 
in a homogeneous, but irregularly shaped volume conductor to be 
evaluated by the summation of surface potentials (Brody, 1957). 
The same principle applies to higher order multipolar equivalents 
also, leading to the general formulation 


anm 2 n—-m ! cos mB 

eae 1 at 4 St « [fren a Bee mgr o) sige 
where g is the specific conductivity of the medium, V, is the sur- 
face potential, and N is the outwardly directed unit normal of the 
surface. Equation 6 is evidently a compact equivalent of that 
which Geselowitz (1960) obtained by a different analytic approach, 
and thus serves to confirm this phase of his study. In the quadri- 
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polar case, equation 6 reduces to 


A209 = g |[vs eas, -ydS, - 2dS,) 
G21 = g |[vs (eas, + ydSx) 


bar = 9 {| Vs(v ds. + 2dSx) (7) 


It 
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boo = 4g Vs (y dS, + 2 dS,) 
bd 


It is now instructive to apply the second degree set of equations 
4 to the case of a dipole of moment, M, located at X, Y, Z. Doing 
so gives 


dao = 2X Mz - YM, -ZM, 

@o1=XM,y+YM, 

bei = XM,+Z2M, (8) 
do2 = > (YM, -ZM,) 

boo = 4 (YM, + Z My) 


Equating the corresponding right sided terms of equations 7 and 8 
gives another set of equations which are essentially identical to 
those which Gabor and Nelson (1954) proposed for locating the 
equivalent dipolar singularity of a source-sink distribution. It is 
apparent from the lead tensor derivation of this particular set of 
Gabor-Nelson equations that their accuracy is vitiated by equiva- 
lent quadripole components, but not by equivalent singularities of 
higher order. Assuming that the equivalent multipolar series is 
convergent, its location can, in principle, be determined by going 
to analogous surface integrals of sufficiently high degree. 


DISCUSSION 


The theoretical developments arising from this study suggest 
that the electrocardiographic contributions of equivalent cardiac 
multipoles can be satisfactorily accounted for by generalizing the 
well-known lead vector principle into a lead tensor apparatus. The 
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study indicates, further, that it is a relatively simple matter to 
cast lead tensor components into a form which relates them directly 
and specifically to the usual zonal and tesseral harmonic notation 
of multipolar singularities. 

Theoretical validation of the lead tensor concept is a fairly 
simple matter if the heart can be contained in a spherical surface 
which does not intersect the body wall. Taking the center of the 
sphere as the coordinate origin, lead fields due to electrodes lo- 
cated external to the spherical surface can be resolved into a 
mathematical series whose components are of the type given in 
equation 4 (Cf. Byerly, 1893, pp, 200-6). Earlier in this report we 
defined the interrelations between actual (rather than equivalent) 
multipolar singularities and the various individual lead field com- 
ponents in terms of lead tensors. According to the principle of 
electrical superposition, the total electrocardiographic effect of 
the multipoles is simply the sum of the individual effects. 

Passing to the case of numerous electrical sources and sinks 
contained within a spherical region, the properties of the individual 
lead field components are inadequate to specify the exact nature of 
the source-sink distribution. Accordingly, the electrical properties 
of the intraspheric milieu can be specified only to the extent of 
their equivalence to current singularities of various order. In order 
to eliminate heart-lead considerations between a multipole of given 
order and lead tensors of higher rank, the equivalent singularities 
are located at the center of the sphere. 

Since usual anatomical relationships admit containment of the 
heart within a sphere, it appears that the lead tensor principle is 
indeed valid, and may be applied to human electrocardiography. 
Implicit in the principle is a confirmation of Geselowitz’s finding 
(1960) that equivalent multipoles can be evaluated in irregularly 
shaped volume conductors. The primary virtue of the lead tensor 
concept, -however, is that it provides a means for completely char- 
acterizing the intrinsic properties of electrocardiographic leads 
with respect to equivalent multipolar representation of the heart. 

The prospects for the accurate experimental determination of 
lead tensors in torso models appear to be excellent. In our labora- 
tory, for example, we have already accumulated considerable ex- 
perience with the precise determination of lead field potentials. As 
we foresee the technical modifications necessary for evaluation of 
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tensor components, our present experimental arrangement would 
probably require little more than the addition of a method for ac- 
curately and systematically transferring the sensing electrode to 
numerous locations throughout the cardiac region of the model. In 
order to facilitate the analysis of many electrocardiographic leads, 
it is likely that the entire technic could be rendered largely auto- 
matic. 

We have found it stimulating to speculate on some of the po- 
tential biomedical applications of the lead tensor apparatus. Cer- 
tainly, a number of interesting possibilities are suggested by the 
concept. Perhaps the most challenging of these is the prospect of 
determining the multipole components of the intact human heart by 
correlative studies involving a human subject and a homogeneous 
torso model patterned from the subject. In principle, this would 
simply be an extension of the method by which Frank (1957) de- 
termined the cardiac dipole components of a human subject. 

Hopefully, this and other highly complex problems will eventually 
yield to future efforts. In the meantime, however, we believe it 
will prove rewarding to pursue less complicated problems, includ- 
ing experimental studies of the lead tensor itself. 


This study was supported by research grant H-1362(C7) of the 
National Institutes of Health, U. S. Public Health Service. 
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A simple model is presented for the formation of functional groups in 
a random neural net. They show the following characteristics: 1. They 
can maintain autonomous activity which might serve as temporary memory 
traces. 2. Early in the process of formation they become resistant to 
contraction. 3. Later they become resistant to expansion. 4. Nearby 
groups inhibit one another. 5. Two groups may contain some cells in 
common. 


There are three assumptions underlying this work. The first is 
that some of the parts of the nervous system that are involved in 
learning processes can meaningfully be treated as though they 
were random nets. The characteristic feature of such a net is that, 
given any two cells in the net, the probability that an axonal proc- 
ess of the first makes a synaptic connection with the second is a 
function only of the distance between the two cells. The mathe- 
matical study of such nets was begun by Shimbel and Rapaport 
(1948), who have since published many papers on the topic. An 
independent treatment has been given by Beurle (1954). Aspects 
of the connective patterns, in particular the dendritic fields, of 
cells in the cerebral cortex, have been discussed by Sholl (1956). 

The second assumption is that learning in the mature mammal 
occurs in two phases. First an experience creates a temporary 
trace within the organism. Then this temporary trace effects a 
permanent change. Properly this should be called a theory rather 
than an assumption, for it has been supported by strong evidence. 
E.g., Duncan (1949) studied the retroactive interference of electro- 
coma on avoidance learning in rats. Some of the control subjects 
were shocked through the legs instead of through the ears; this 
type of shock seemed to be more traumatic but it did not produce 
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unconsciousness and it caused relatively little interference. Dun- 
can’s data indicated that it took about fifteen minutes for the per- 
manent effect of a single experience to be fully established. 

The third assumption concerns the nature of the permanent change. 
There does not seem to be any good evidence on this matter. A 
variety of hypotheses have been proposed. We shall adopt the 
suggestion made by Hebb (1949), Eccles (1953) and others, that 
the locus of the permanent change is at the synapse: if one cell 
repeatedly fires just before another, with which it has a connec- 
tion, then the influence of the first upon the second becomes 
more excitatory. 

If the permanent change is of this nature, then the temporary 
trace must consist, at least in part, of persistent activity in par- 
ticular loci. The mechanism that has often been proposed (e.g., 
Rashevsky, 1937) for such a temporary trace is the self re-exciting 
cycle: for example, cell A excites B, B excites C and then C ex- 
cites A again. This example is over-simplified; in general, a cell 
needs to be stimulated by a number of others within a sufficiently 
short interval of time in order to fire. Instead of a simple chain, 
a group of cells with many interconnections is required. It has 
been shown in isolated slabs of the cerebral cortex of the cat that 
some cells can maintain activity for short periods by this mecha- 
nism of mutual excitation (Burns, 1958, page 28). An aggregate of 
cells capable of maintaining activity has been given the name 
‘cell assembly’’ by Hebb (op. cit.). Hebb uses this as a central 
concept in a broad scheme of psychological theory. He cites evi- 
dence that in the sensorily deprived mammal (that presumably has 
not been able to develop the needed cell assemblies) learning is 
very inefficient. 

There are certain difficulties with the concept of the cell as- 
sembly. Beurle (op. cit.) has shown that activity will tend to 
spread in waves across a random net, because of the refractory 
period of the cells. Rosenblatt (1958) argues that if there were a 
number of cell assemblies within a net, with the passage of time 
these assemblies would coalesce, by Hebb’s hypothesis of in- 
creasing synaptic efficiency. Milner (1957) has suggested a 
scheme for the development of cell assemblies that meets these 
objections. It is the purpose of the present work to explore mathe- 
matically a simpler model. This model, in common with Milner’s, 
involves inhibitory as well as excitatory processes. Evidence on 
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inhibition in the brain stem and cerebral cortex has been given by 
Purpura (1958). 

Let x be a vector denoting position in a random net. Let f(z) 
be the input at z: the rate at which impulses from cells outside the 
net are impinging upon a cell at 2 Let g(x) be the output at 2: 
the rate at which a cell at « is firing. Let (zx) be the net excita- 
tion at z, due to both the input and the activity of cells within the 
net. Each of these quantities, f(z), g(x) and €(z), is a local aver- 
age, each varies with time and each could have the units of im- 
pulses per second. Let p(z) be the density of cells at 2. Finally, 
let @(z,zx’) be a parameter denoting the influence of activity at 2’ 
upon acell at z. Then 


e(2) = F(a) + | o(e") g(a’) 0(@, 2°) de’ (1) 
where the integral extends over the net. 

In general 6(z,2’) depends upon two factors. The first is the 
distance between z and x’. In this paper, however, we shall be 
concerned with a few groups of cells in neighboring or overlapping 
regions and the spatial factor will not be of importance. The other 
factor, which will be important here, is the amount of nearly simul- 
taneous activity that has taken place at 2 and x’. Following such 
activity, @(z,z2’) increases by Hebb’s hypothesis. 9 is a mono- 
tonic function of time. If there has been no simultaneous activity 
at 2 and x’, @(z,2’) has its least possible value. Since neural 
inhibition is known to be important, we suppose this minimum value 
to be negative. The idea that 6 can change from a negative to a 
positive value does not imply that any single cell changes the 
nature of its influence from inhibitory to excitatory. It means that 
initially the net effect of activity of the one region upon the other 
is inhibitory and that with simultaneous activity the excitatory con- 
nections can be permanently potentiated so that this net effect be- 
comes excitatory. The initial, minimum value of @ will be denoted 
by -j. It does not seem reasonable for @ to increase without 
bounds; we suppose @ to have a maximum attainable value, Om. 

Essentially what happens in this model is that the types of ac- 
tivity stimulated by various inputs become increasingly stereo- 
typed and fixed. As will be shown below, if a particular input 
pattern has been prolonged and repeated enough to cause sufficient 
‘training’? among the stimulated cells, then any similar input, 
even though it may stimulate more or fewer cells, will cause the 
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same pattern of activity; the cells that have been ‘‘trained’’ to- 
gether (that is, have previously had persistent simultaneous ac- 
tivity) will now stimulate one another but inhibit other cells in the 
vicinity. An input pattern that is continuously graded over space 
tends to elicit activity with sharp demarcations. We shall discuss 
discrete groups of cells; our results will show to what extent the 
boundaries of the groups are actually fixed and to what extent they 
are-variable. If the cells of the zth group have been trained to- 
gether, but none have been trained with any cells in the other 
groups being considered, equation (1) can be rewritten, 


Gi = fit Oinigi-j) mv gv, (2) 


vei 


where n; is the number of cells in the th group. 

It is conceivable that the variation in time of f(z) might be too 
sharp to permit the development of anything like a steady state. 
However, the analysis of the present paper will be limited to such 
relatively simple conditions. In the steady state, the output at any 
place is a certain function of the net excitation there: g = ¢(€). 
This function has been discussed by Rapaport (1950). It depends 
upon the following properties of the individual cell: the refractory 
period 5, the period of temporal summation oa, and the threshold A. 
If 4 2 2, the curve has a sigmoid flexure; for large values of ¢€, the 
curve is asymptotic to 1/5. The inverse function é 1'(g) gives 
the amount of excitation needed to maintain the output at a level! 
g- In the discussion below we shall make use of the following 
linear approximation: 


6x10? 
4x103 


2x10 


2" (impulses per second) 


0 25 50 75 100 
g (impulses per second) 


FiaurE 1. ¢~'(g), with the linear approximation sg +k, for th. 
= 1072 sec., = i073 sec., andAh=5. g » for the case 


CELL ASSEMBLIES 47 
~ *(9)=3g9+k, (3) 


where 0<g<1/6. In order to get a rough idea of the dependence 
of s and & upon the parameters 5, o and h, & '(g) was plotted for 
various values of these parameters and it was found that s is 
about 2 (28/c)°-7? and & is h?*3/408° 399-7, 

The first case to be discussed will be the simplest: a group of 
nm cells which have been trained together, and none others in the 
vicinity, are stimulated by an input f. Now 


E=f+Ong. (4) 


In Figure 2a, f + Ong is superimposed upon ¢ !(g). The two 
curves intersect when g=go. If g< go, then f+ Ong>¢ '(9); 
that is, the excitation € is greater than the excitation needed to 
maintain the output g and g will increase. If g > go, then f+ Ong < 
@ ‘(g); that is, the excitation € is less than the excitation needed 
to maintain the output g and g will decrease. So there is stable 
equilibrium at go. 

When the input ceases, that is, when f falls to zero, then g, 
falls to zero. But if the input is sufficiently prolonged and re- 
peated, the parameter 9 may attain a value great enough that there 
are two more points of intersection of Ong and ¢ ‘(g). It can be 
seen in Figure 2b that zero and g’”’ are values of stable equilib- 
rium and that g’ is a threshold. That is, if g < g’ then the output 
drops to zero, but if g > g’ then the output moves toward g’’. Sup- 
pose the group to be initially at rest, with g = 0; an input f can 
push g over to its other stable value g’’ and even after f falls back 
to zero the group will maintain this level of activity, until it is af- 
fected by inhibition or fatigue. Rapaport named this the ‘‘ignition”’ 


f+60 


Jo 
FIGURE 2a. FIGURE 2b. 
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phenomenon of random nets (1952). The product @n of the present 
discussion corresponds to the parameter a, the axon density, in 
Rapaport’s paper. The least value of @n such that @ng has a non- 
zero intersection with ¢ 1(g) is the least value such that ignition 
is possible. This critical value will be denoted by c. The prob- 
lem of determining c as a function of 4 was studied by Trucco 
(1952) for a model in which the cells could act only at unit inter- 
vals and there was no temporal summation. In the present model, 
c~s+ké. 

The sizes of the assemblies that develop depend upon the in- 
puts f. We might try identifying f with, for example, a sensory 
after-discharge, but there does not seem to be pertinent evidence. 
™ this model, the ease with which an assembly can be formed 
varies with the size of the group. On the one hand, the smaller 
the number n, the greater 6 must become in order that 6n2 c, the 
condition that the group be ignitable. (The least value of n such 
that this condition can ever be met is c/@m. For the case illus- 
trated in Figure 1, c = 36; if for example, 6m = 0.1, then the small- 
est possible assembly would contain 360 of these cells.) On the 
other hand, the greater the number n, the more the initial inhibi- 
tion, —jn, that f has to overcome in training the assembly. 

In the second case to be considered, n; + ng cells have been 
given some degree of training together and only n, of the cells are 
stimulated by an input f. At equilibrium, 


& '(91) =f + Oni 91 + Onaga, (5) 


b *(92) = On1 91 + One go. (6) 
When the input f drops to zero, the excitation is the same through- 
out. Therefore, as long as @ is uniform, the n, cells cannot be ig- 
nited unless the whole set of n; + no is. Only if there has been 
so little training that @ is still negative and activity of the ny, 
cells inhibits the nog cells can an assembly form of the n; cells 
only. 
In the third case n; cells are trained together and an input f 


stimulates mi +g cells simultaneously. Using the linear ap- 
proximation, the equilibrium conditions are 


Sgit+k=f+ 010191 - jnega, (7) 


SJ2+k=f-—jnigi —jnage. (8) 
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The solution of this system of equation is 


91 = (f-&)s/A, 
G2 =(f-k)[s- (61 + fynil/A, (9) 
A = (8 - 8111)(8 + jna) - 7? m1 no. 


According to this solution, g2 < 0 if 6101 >s — jn; (assuming that 
f-— and A are positive). Of course g2 cannot be negative. This 
means that if 611 > s — jni, the ng cells are inhibited and g2 = 0. 
Thus if there has been sufficient training, the assembly is resist- 
ant to expansion. However, just as the linear approximation is 
not valid when it gives gg a negative value, neither is it valid 


when it gives g; a value greater than 1/5. If f is very large, g1 
approaches 1/5; then, from eauation (8), 


eet 175) 

f2 S+JNg 4 (10) 
Thus no matter how well trained the n; cells are, no matter how 
large 61 is, a sufficiently large f can force simultaneous activity 
in the no cells. If this simultaneous activity is sufficiently pro- 
longed and repeated, it will cause incorporation of the ne cells 
into the assembly. 

To summarize what was found in the second and third cases: 
After a group of cells has had enough activity that @ > 0, it is im- 
possible for an assembly to be formed from only a part of the group. 
Later in the training, when @ attains the value (s/n) - j, it be- 
comes difficult, although not impossible for a larger assembly to 
form with the given group as a part. 

In the fourth case two groups of cells, of sizes n; and nz, have 
been trained separately and receive an input / simultaneously. 
The approximate equilibrium conditions are 


89, +k=f+ O01inigi1—jnege2, (11) 
8Jo+k=f—jnigi+t 027292. (12) 

Solving (11) for g1 and (12) for ga, 
Ge (13) 


8s-—O,n, 


a ae 
Geol k-jm 91 


8S — O02n2 


(14) 
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(Of course these equations have no utility if 62 s/n for either 
group.) Considering g1 in equation (13) to be a function of ga, it 
can be seen that the value of 6; affects the gi-intercept and the 
slope of the function, but that (as long as 41 < s/n1) the ga-inter- 
cept is unaffected: if go 2(f—)/jna, the m1 cells receive at 
least as much inhibition as excitation and g; = 0. Similarly, when 
ga in equation (14) is considered to be a function of gi, the value 
of 92 affects the go-intercept but not the g,-intercept. All of this 
is pictured in Figure 3. In 3a, 61 and 62 are both small; in 3b 
they have larger values. 

Suppose gi and gz to have values such that the point (g1,92) 
lies above the graph of equation (13). This means that the excita- 
tion of the n; cells (the right side of equation (11)) is less than 
that required to maintain g; (the left side of equation (11)) and g1 
will decrease. If (gi, 92) is below the graph of (13) there is more 
excitation than needed and g; will increase. A similar analysis 
applies to g2. This is illustrated by the arrows in Figure 3. Only 
the circled points are points of stable equilibrium. The condition 
that there be a stable point where neither g; nor g2 vanishes is 
that 6< (s/n) -j for both groups. If that is the case, prolonged 
simultaneous stimulation can cause consolidation of the two groups. 
But if the groups have had enough training, activity in one can 
completely inhibit the other. The group with the larger value of On 
tends to monopolize the activity. Alternatively, if @n is about the 


"2 Be alk) 
f-k 
fh 1 
ns ila) 
9, 
fk. fok 
S-6,N, n, 


FIGURE 8a. 
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FIGURE 3b. 


same for the two groups but f is not uniform, the two groups form a 
mechanism somewhat like that proposed by Landahl (1938) for 
discrimination. 

The situation is much the same if the two groups are ignitable. 


Let 61 = 02 = Om, 21> ng and Om M2. If gi and g2 were almost 
1/5, then 


Sete avg Pana gi fan a Ih (15) 
6 


The activity of a group can remain at almost 1/5 only so long as 
the excitation is not less than(s/5)+k=c/6. If inhibition causes 
the excitation to fall below this, the activity decreases, the exci- 
tation drops further and if the inhibition continues the activity 
vanishes. So activity can continue in the second group when f = 0 
only if 


(16) 


If assemblies are to endure without too much coalescence, the 
above condition should be unlikely even when the two groups are 
very large and about the same size. This implies that 0m Sj (t.e., 
that the final value of @ is not greater than the absolute magnitude 
of its initial negative value.) 

The discussions of cases two, three and four imply that the cells 
of a given assembly act all together or not at all, that the assem- 
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blies are atomic. The fifth case shows an exception to this. Sup- 
pose that 1 + m2 cells have been trained together and that na +13 
cells are stimulated simultaneously with an input f. That is, the 
nz cells constitute the intersection of the set of trained cells with 
the set of cells receiving the input. By equation (2) the values of 
the excitation are as follows: 

E, = O(n1gi + n2g2)—-jn3 93, 

‘Eg =f + 6(n1g1+ 292) -—j73 93, (17) 

Eg =f—j(migitm2g2+ nz). 
If Ono > jng, then @n2g2 > jnsga; the ni cells will be activated 
and if f is not too large the ng cells will be inhibited. On the 
other hand, the n; cells may receive more inhibition than exci- 
tation. To investigate this possibility, we set g1 = 0, 2 = sg2 + 
k and €3 = sg3 + k, solve for g2 and g3 (as for g: and g2 in case 
three) and set On2 g2 < jn3g3: 


Ono (f — k)s/A < jna(f— k) ls — (6 + j)nal/A. (18) 
Solving the inequality for no, 
-1 
6 64+ 7 
KEYS || See i < 
2 E — | (19) 


If (19) is satisfied, then upon the application of this type of input 
the n; cells will be inhibited and eventually an assembly can be 
formed of the no + ng cells. Let the training parameter of this as- 
sembly be 0’. If @’n, < jn, it is still possible for the n; + M2 as- 
sembly to be activated without activating the no + ns assembly. 
Thus two assemblies can contain a limited number of cells in com- 
mon without losing their separate identities. 

In the opinion of the author, an attempt directly to relate the 
phenomena discussed in this paper with psychological evidence 
would be misleading. If there are such things as cell assemblies, 
evidence relating to the ‘‘law of mass action’? indicates that it is 
the interaction of great numbers of them that influences behavior. 
There are about 10'° neurons in the nervous system; if one in ten 


were incorporated into an assembly and the mean assembly size 
were 104, there would be 105 assemblies. 
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Quantum biology is the quantum mechanical study of electrons in 
molecules of biological interest. This requires the solution of problems 
involving many electrons. Approximation methods are therefore neces- 
sary and are discussed. The present study, concerned with the mobile 
electrons in riboflavin (FMN) and its radicals (FMN , FMNH and 
FMNH,+), is based on the approximation method, developed by B. Pullman 
and A. Pullman. The solution of the eigenvalue problems so obtained 
gives the energy levels of the mobile electron systems involved. The 
corresponding eigenvectors yield the mobile electronic charges of the 
atoms of riboflavin radicals which have contributed mobile electrons. 
Important differences of the net charge distributions of these radicals 
are emphasized. The longest wave length of light absorption is calcu- 
lated from the obtained energy levels and agrees, within the accuracy of 
the method, with corresponding experimental results. From the appropri- 
ate calculated results, electronic assignments are obtained for the ex- 
perimental transitions involved. 


Introduction. The fundamental process on which all life on earth 
depends is the absorption of sunlight by pigments like chlorophyll 
of green plants. By this absorption, a mobile electron of the 
chlorophyll molecule is raised from one discrete energy level to 
another. 

In the respiratory oxidation-reduction chain, mobile electrons 
are transferred from a discrete energy level of diphosphopyridine 
nucleotide to another in riboflavin, and then to other levels in the 
cytochromes. The study of the electrons in these molecules on 
the basis of quantum mechanics is of fundamental biological im- 


portance, 
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In this paper, a study will be presented of the mobile electron 
systems of riboflavin and of the radicals which arise from the 
latter by the acceptance of an electron and one or two protons. 
As will be seen below, properties of the individual atoms and 
bonds of the mobile electron system can be taken into account ex- 
plicitly. Such studies of electrons in molecules of biological in- 
terest, based on quantum mechanics and taking into account ex- 
plicitly the individual atoms and bonds of the molecules, may be 
said to belong to the field of quantum biology. Thus defined, 
quantum biology is a part of submolecular biology. It should be 
distinguished from other applications of quantum mechanics to 
biology, e.g. target theory. 

In many molecules of biological interest and in all the above- 
mentioned molecules there are alternate single and double bonds 
assigned by classical chemistry which are known as conjugate 
systems. Quantum chemistry (B. Pullman and A. Pullman, 1952; 
C. A. Coulson, 1952; K. S. Pitzer, 1953; W. Kauzmann, 1957) 
shows that this alternation is, to a large extent, only apparent be- 
cause the mobile electrons are actually non-localized. The ease 
of redistribution (i.e., delocalization) of these mobile electrons 
makes their behavior important for long range effects, such as oc- 
cur in biology. Since these mobile electrons in general are much 
less tightly bound than the localized electrons, their excitation 
energies are smaller than those of the latter. Consequently, the 
excitation of the former accounts for the longest wave absorption 
peaks of such systems. It is known from x-ray studies that the 
conjugated systems are planar, This is explained in quantum 
chemistry (Pullman and Pullman, op. cit.; Coulson, op. cit.) from 
the principle of maximum overlap (L. Pauling, 1948). 

As is well known, the behavior of an electron either within or 
outside matter is described in theoretical physics by the Schré- 
dinger equation, Solution of this equation involves an eigenvalue 
problem determined by a partial differential equation with appropri- 
ate boundary conditions. In general, it is not yet practical to ob- 
tain rigorous solutions of the Schrédinger equation for problems of 
more than one electron because of the interactions between the 
electrons, 


Therefore, to make quantitative progress it has been necessary 
to apply approximation methods such as the Valence Bond (V.B.) 
the Molecular Orbital (M.O.), or the Free Electron (F.E.) meted 


QUANTUM BIOLOGY 57 


(Pullman and Pullman, op. cit.; Coulson, op. cit.; Kauzman, op. 
cit.). Moreover, in the treatment of heterocyclic systems, certain 
experimentally determined quantities must be known (R. Pariser 
and R. G. Parr, 19538). 

Pullman and Pullman (op. cit.) have developed a method for the 
approximate wave mechanical treatment of a system of mobile 
electrons, which has led to many fruitful results in the study of 
carcinogenesis (A. Pullman and B. Pullman, 1955), of biology in 
general (for a recent review, B. Pullman and A. Pullman, 1960), 
and more recently of neuropharmacology (G. Karreman, I. Isenberg 
and A. Szent-Gyorgyi, 1959). This method is based on the Linear 
Combination of the Atomic Orbital approximation of the Molecular 
Orbital Method. 

In the latter approximation, a molecular orbital wave function, 
ww, is developed in a linear combination of the mobile electronic 
atomic orbitals 9, of all (n) atoms participating in the conjugation: 


n 


=) 6,9; (i =1---, 7) (1) 


t= 1 


The coefficients c, form the eigenvector ¢. They are obtained as 
the solution for the eigenvector corresponding to the eigenvalue 
which gives the energy E of the M.O. The method of determination 
of these energy levels and the corresponding eigenvectors will be 
outlined now. 

Application of the Ragleigh-Ritz variational method leads to the 
equations (Pullman and Pullman, 1952; Kauzman, 1957); 


>, 4, - £4) c,= 0(j =1---n) (2) 
i=! 
Here the integrals Hj, called Coulomb integral when ¢ = j and ex- 
change integral when z # j, are given by: 


H;, - | 949, dv (3) 

H being the Hamiltonian operator and the overlap integral S, by: 
s 4 

S,. -/ 9, ,de (4) 


both integrated over all the coordinates involved. The Coulomb 
integral H,, is a measure of the electro-negativity of atom @ for a 
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mobile electron and is expressed in the method of B. Pullman and 
A. Pullman (1952, 1958a) as 4 +6,8. The parameters a and B 
are H,, (each atom 7 being a carbon atom) and H,,, in which the 
atoms i and j are both a carbon atom, respectively. Both of these 
parameters have negative values (Coulson, op. cit.). In the same 
method, the exchange integrals A; are expressed in terms of 8 as 
unit: H,, = p; ;B- The values of the parameters 6, and Pi j (Pullman 
and Pullman, 1958a) are, therefore, measures of the electro- 
negativity of atom ¢ for a mobile electron, and of the exchange of a 
mobile electron between atoms 7 and j respectively. They are ob- 
tained from experimental data such as dipole moments, bond 
lengths and nuclear electric quadrupole coupling parameters (Pull- 
man and Pullman, 1952; J. H. Goldstein, 1956). 

The system of equations (2) has a non-trivial solution only when 
the secular equation (5) is satisfied. 


|H,, -8,,E | =0 (5) 


In the approximation method of Pullman and Pullman (1958a), 
the following simplifications are made: S, =1l if z= 7, and 0 if 
tet 4, H,, is only #0 if ¢ and j refer to nearest neighbors. The 
latter expresses that only the exchange of mobile electrons be- 
tween nearest neighbors is taken into account. The solution of 
this equation for E in the case of a molecule, in which n atoms 
participate in the conjugation, yields the energy eigenvalues 
E.(j=1---n) in the form BE, =%+ K,8 (Pullman and Pullman, 
1952; 1958b). Therefore, % is important only in determining the 
absolute level of the energies E., while the relative values of E. 

; j j 
are determined by the value of 8. Substitution of such an elgen- 
value E; into equations (2) leads to n linear compatible equations 
for the n unknowns c.... Their solutions are the components of the 
elgenvector c;, corresponding to the eigenvalue E.. From the latter 
the mobile electronic charge e, of the ¢*® atom participating in the 
conjugation (contributing electrons to the mobile electron system) 
is given by (Coulson, op. cit.): 

n 
e, =>° a,(c; .)? (6) 
j=1 
Similarly, the mobile order (or mobile electron strength), P,, of the 
bond between atoms i and &, both participating in the conjugation, 
is obtained from (Coulson, op. cit.): 
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n 


Pin = 215 %5%,; (7) 


j=l 


In the latter formulae, a; represents the number of electrons which 
occupy the M.O. w,- Therefore, a; = 0,1 or 2 for ¥,, being occupied 
by zero, one or two electrons, respectively. 

Though this method contains several simplifications as outlined 
above, I. Isenberg and A. Szent-Gyorgyi (1959) have shown that it 
can be applied to the calculation of spectra of several molecules 
of biological interest, including riboflavin (or FMN*). This method 
is used in the present work. It extends the earlier work of G. 
Karreman on riboflavin, first presented elsewhere by Pullman and 
Pullman (1958b, 1959). In this last paper Pullman and Pullman 
published their important finding of the anti-bonding character of 
the highest filled mobile electronic molecular orbital in FMNH,. 
A similar result has recently been found for the tranquilizer 
chlorpromazine with obvious implications for neuropharmacology 
(G. Karreman, I. Isenberg and A. Szent-Gyorgyi, op. cit.). 

Some General Considerations of the FMN-Semiquinones System. 
It is known (Pullman and Pullman, 1958b) that FMN is a good 
electron acceptor. We will consider a molecule of FMN which has 
accepted an electron, e.g. from a reducing agent, in solutions of 
varying pH. Besides FMN, the semiquinones FMN~, FMNH and 
FMNH; will be present in such solutions depending on the pH. 
These species participate in the following equilibria: 


FMN +e <2 FMN- (8) 
in neutral solution, 
FMN + e~ +H* = FMNH (9) 
in weak acid solution, 
FMN +e +2H* = FMNHZ (10) 


in strong acid solution. 52 ge 
The constants of the last three equilibria will be called Kk, K, 
K. respectively. Furthermore, the relation 


and K, 
[FMN] + [FMN7] + [FMNH] + [FMNH2] =[FMN], (11) 


*Since riboflavin and flavin mononucleotide, FMN, have the same sys- 
tem of mobile electrons, they will both be designated by the abbreviation 
FMN in this paper. 
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in which [FMN], is the total riboflavin concentration, has to be 
taken into account. ‘Then the following relations can easily be 


derived in the case of dilute solutions (in which the mite, a 
FMN~ 


approximately equal to the concentrations) for ” ~ (EMNIS 


[FMNH] [FMNH#] 
eNO) 
Y [FMN], [FMN], 


erie re (12) 
1 +(a/p) + ax + aga” 


ax 
ae ee ee (13) 
1 +(a/p) + ax + aqz 


2 
yl Bae cvkies bok oe (14) 
1 + (a/p) + ae + aga? 


In these formulae « =[H"], a = UTA Rog ie g = i. so that z, 
a, p and qg are always positive. Therefore it is easy to see that 
with increasing [H*], [FMN7] will steadily decrease; [FMNH] will 
go through a maximum and [FMNH3'] will steadily increase. 

Experimentally, it is known (H. Beinert, 1956b) that a neutral so- 
lution containing only riboflavin has its longest wave absorption 
peak at 445 mu. However, in the presence of a reducing agent, the 
longest wave absorption peak of riboflavin is at 565 mu in weak 
acid solution and at 503 mu in strong acid solution (H. Beinert, op. 
cit; A. Szent-Gyorgyi, 1960). 

Mobile Electronic Structures of FMN and its Semiquinones. The 
mobile electron system of FMN is that of isoalloxazine, represented 
in Figure 1, together with the numbering of its atoms participating 
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FIGURE 1, 
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in the conjugation, H, by hyper-conjugation (Pitzer, 1953). This 
System contains 22 mobile electrons, obtained as follows. Each 
of the C and O atoms, as well as the N atoms number 1 and number 
5, contributes one electron to the mobile electron system. How- 
ever, N atoms number 3 and number 10 and each of the hyper- 
conjugated CH; groups contribute two electrons to the mobile 
electron system. 

Therefore, taking into account the Pauli principle, each of the 
11 mobile electronic molecular orbitals with lowest energy levels 
is occupied by two electrons with opposed spins, in the ground- 
state of FMN. For FMN the highest calculated energy level of the 
filled and the lowest energy level of the empty mobile electronic 
molecular orbitals are given in Figure 2. The energy levels are 
counted from the lowest of the energy levels of the mobile elec- 
tronic molecular orbitals. The energy levels were obtained by 
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FIGURE 2. Highest (= 11th) energy level of the filled and lowest 
(= 12th) energy level of the empty mobile electronic molecular orbitals 
of FMN. 


solving an eigenvalue problem involving a matrix of the order 20 by 
20*. Electrons are represented in the figures of this paper by 
small arrows indicating their spins. In Figure 3 are reproduced 
the calculated net contributions of the delocalized mobile elec- 
trons to the electrical charges of the atoms in FMN (the electronic 
charge being the unit). Such a contribution, called the net or formal 
charge, of an atom participating in the conjugation, is obtained in 
the following way. The mobile electronic charge of that atom is 
subtracted from the number of electrons contributed by that atom to 
the mobile system. These net contributions of the mobile electrons 
to the charge of an atom will be given in this paper. The mobile 


*As is nearly always the case for molecules of biological interest, the 
solution can be obtained only by the use of a large digital computer, for 
_ which a Univac is used in this work. 
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electronic charges of the atoms in FMN are given in Figure 7 of 
the paper by Pullman and Pullman (1959). In the figures of the 
present paper which give the net charges of the atoms, only the 
relative positions of the atoms and no bonds are given. The 
reason for such a rather unusual representation is that it is in- 
consistent to use simultaneously both electric charges and classi- 
cal bonds side by side in molecular diagrams (K. Fajans, 1949; 
1959). 

The l(owest) e(nergy) l(evel) (of the) e(mpty) (mobile electronic) 
m(olecular) o(rbital) determines the molecular electron affinity. 
Comparison of the value of the l. e. 1. e. m.o.’s of molecules of 
biological importance shows that riboflavin is a good electron ac- 
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FIGURE 3. Net contributions ofthe delocalized mobile electrons to the 
net electrical charges of the atoms in FMN. 


ceptor (Pullman and Pullman, 1958b, 1959). After accepting an 
electron, FMN becomes FMN™. In the ground state of the latter, 
the accepted electron occupies the 12th molecular orbital which is 
the empty orbital with the lowest energy level in FMN. In Figure 
4 are represented the highest and lowest energy levels of the filled 
and empty mobile electronic molecular orbitals respectively in 
FMN’. In Figure 5 are given the net mobile electronic charges of 
the atoms in FMN , resulting from delocalization of the 23 electrons 
in its mobile system. 

From Figure 5 it is seen that of all the atoms in FMN , the N 
atom number 1 has the greatest negative net charge. Consequently, 
the H* ion which will be taken up first by FMN’, that is the H* ion 
involved in the equilibrium (9), will be taken up by N atom number 
1. As a consequence, the valence state of N atom number 1 is 
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FIGURE 4. Highest (= 11th and 12th) energy levels of the filled and 
lowest (= 12th and 13th) energy levels of the empty mobile electronic 
molecular orbitals of FMN~. The 12th mobile electronic molecular orbital 
is half filled and half empty. 


changed, as the two lone pair of electrons of the latter atom now 
participate in the conjugation. Because of this finding, the radical 
FMNH (Figure 6) with 23 electrons in its mobile system is also 
studied. It is emphasized here that the radical illustrated in 
Figure 6 is more stable than the one which arises if the H® ion had 
been taken up at N atom number 5 (the latter is probably not even 
planar). The preferential fixation of the H* ion in the 1- position 
of FMN and FMN’ is consistent with all the Ultra Violet spectro- 
photometric work of P. Hemmerich (1960, personal communication) 
of the University of Basel, Switzerland. 

The radical FMNH has 23 electrons in its mobile system. There- 
fore in its groundstate the molecular orbitals corresponding to the 
lowest 11 energy levels will be occupied by two electrons, and the 
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FIGURE 5. Net contributions of the delocalized mobile electrons to the 
net electrical charges of the atoms in FMN . 
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FIGURE 6. 


one corresponding to the next higher level by one electron. In 
Figure 7 are given the results of the calculations for the 11th, 
12th and 18th levels. From the & value of the 12th energy level it 
is seen that FMNH is an excellent electron donor. It is also a 
good electron acceptor. Because of the simultaneous occurrence 
of these two properties this radical FMNH may be a very important 
intermediate in the biological reduction and oxidation of FMN (H. 
Beinert, 1956b). Because of these properties it can probably form 
a self charge transfer complex (FMNH)., as pointed out by P. 
Hemmerich (1960, personal communication). In Figure 8 are repro- 
duced the calculated net charges of the atoms participating in the 
conjugation of FMNH. It is important for later reference to notice 
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FIGURE 7. Highest (= 11th and 12th) energy levels of the filled and 
lowest (= 12th and 13th) energy levels of the empty mobile electronic 
molecular orbitals of FMNH. The 12th mobile electronic molecular or- 
bital is half filled and half empty. 
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FIGURE 8. Net contributions of the delocalized mobile electrons to the 
net electrical charges of the atoms in FMNH. 


that in FMNH, N atom number 1 has a great positive net charge as 
opposed to the great negative net charge of that atom in FMN . 
Another proton can be taken up by FMNH resulting in the radical 
FMNH3, which has the same arrangement of its atoms as in the 
molecule FMNHg2 (see the corresponding Figure on page 137 in the 
paper by Pullman and Pullman, 1959). However, the radical 
FMNHz2 has only 23 electrons in its mobile electronic system com- 
pared with 24 in the case of the molecule FMNHe. Therefore, in 
the groundstate of this radical the occupation numbers of the 
levels are the same as those of FMN’ and FMNH. In Figure 9 are 
represented the calculated 11th, 12th and 18th energy levels of the 
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FIGURE 9. Highest (= 11th and 12th) energy levels of the filled and 
lowest (12th and 18th) energy levels of the empty mobile electronic 
molecular orbitals of FMNH¢+. The 12th mobile electronic molecular 


orbital is half filled and half empty. 


energy 
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mobile electron system of FMNH;. In Figure 10 the net charge dis- 
tribution is reproduced as calculated for the radical FMNB;. Again 
for reference in later work it is noted that in FMNH3 also N atom 
number 5 has a great positive net charge as opposed to the great 
negative net charge of that atom in FMN and FMNH. 

From the difference of any two energy levels the corresponding 
wave length of the transition can be calculated using the value 
B = 8.26 e V, obtained by I. Isenberg and A. Szent-Gyorgyi (1959). 
In Figures 2, 4, 7 and 9 the values so obtained for the wave 
lengths in mp are represented along the dotted arrows. From these 
values, the following electronic assignments are obtained for the 
absorption peaks of the semiquinones of FMN in weak and strong 
acid solutions. 

The lowest energy transition of FMN in neutral solution has an 
experimentally observed wave length of 450 mu. FMN is one of the 
substances used by Isenberg and Szent-Gyorgyi (op. czt.) to calcu- 
late the value 3.26 eV for 8 and the results for riboflavin in par- 
ticular correspond closely to this average value of B as can be 
seen from their table (1). Consequently, the longest wave length 
of the absorption peak of FMN calculated with this value of B 
corresponds closely to the observed wave length. The experi- 
mentally observed longest wave lengths of the absorption peaks of 
the semi-quinones of FMN are quoted by H. Beinert (1956b), 
565 mp in weak acid and 503 mp in strong acid solution. These 
are to be compared with the following theoretically calculated 
wave lengths of the lowest energy transitions (which, of course, 
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FIGURE 10. Net contributions of the delocalized mobile electrons to 
the net electrical charges of the atoms in FMNH¢. 
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correspond to the longest wave lengths of absorption peaks): 
FMN 565 my; FMNH 576 mp; and FMNH+ 480 mp. From this com- 
parison and taking into account the accuracy of the method (I. 
Isenberg and A. Szent-Gyorgyi, 1959), the following electronic 
assignments are proposed. (1) The absorption peak at 503 mp is 
due to the transition of an electron from the 11th to 12th mobile 
electronic molecular orbital in FMNH2. (2) The observed wave 
length of maximum absorption at 565 my in weak acid solution is 
probably due to the transition of an electron in FMNH from the 12th 
to the 13th mobile electronic molecular orbital. Although the long- 
est wave length of maximum absorption in FMN” is at 565 mp, the 
latter assignment is made because FMNH is the dominating species 
in weak acid medium. Nevertheless, the observed peak at 565 mp 
may be due to a mixture of FMNH and FMN . 

In conclusion, this work shows that (1) the net charge distribu- 
tions of FMN , FMNH and FMNH¢ reveal the following important 
features. In FMN the atoms number 1 and 5 have a great nega- 
tive charge. In FMNH the N atom number 1 has a great positive 
charge and N atom number 5 a great negative charge. In FMNH+ 
the N atoms number 1 and 5 have great positive charges. (2) The 
calculated energy levels of FMN , FMNH and FMNH¢? allow the 
assignment of (a) the experimentally observed absorption peak at 
503 mp to the transition of an electron from the 11th to the 12th 
mobile electronic molecular orbital in FMNH¢ and (b) the experi- 
mentally observed absorption peak at 565 mp to the transition of an 
electron from the 12th to the 13th mobile electronic molecular 
orbital in FMNH (and in FMN ). 
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A postulated scheme for sequential enzyme synthesis coupled with 
linear chain reaction is proposed. Mathematical calculations are carried 
out for the transient-state formation of enzyme coupled with enzyme- 
substrate reactions. Examples are presented to demonstrate the relations 
between various operational parameters during enzyme synthesis as a 
function of time. Significance of the phenomenon of sequential synthesis 
is discussed in relation to abnormal growth and cellular injury. 


In previous publications (Heinmets and Herschman, 1960; Hein- 
mets, 1960) models for the induced enzyme synthesis were proposed 
and mathematical analysis was carried out. It was pointed out that 
in spite of the fact that the role of the inducer is not quite certain, 
there is sufficient experimental evidence to show that sequential 
induction systems are operative in cells (Stanier, 1947 and 1955; 
Suda ez al; 1950a,b) and such systems lend themselves for an anal- 
ysis. The subject of sequential enzyme induction has been thor- 
oughly reviewed and analyzed in the literature (Pollock, 1959; Dixon 
and Webb, 1958). At present we shall attempt to carry out a formu- 
lation and mathematical analysis of sequential enzyme formation in 
which a derepressor role for the inducer is indicated. We are aware 
of limitations of such work and perhaps simplicity of model-system 
may seem unrealistic. However, when we proceed further in our 
work, it will become apparent that mathematical treatment of even 
such a simple system becomes quite complex and a more complete 
analysis is warranted only when systems have been experimentally 
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elucidated so that a more definite model-system can be proposed. 
Nevertheless, we hope that even an analysis of such a simple sys- 
tem may induce some quantitative experiments in sequential enzyme 
synthesis. 

In previous publications (Heinmets and Herschman, 1960; Hein- 
mets, 1960) a number of postulates have been made concerning the 
basic operational units of the model system. These are not re- 
peated here and these publications should be consulted. 


Sequential Synthesis of Enzymes 


A postulated relation* between enzymes and substrates in chain 
reactions as well as in synthetic processes is presented schemati- 
cally in Fig. 1. Here the role of substrate is twofold. First, it 
interacts with an enzyme in a transfer reaction and, secondly, it 
acts as an inducer in the synthetic system for substrate specific 
template activation. In order to carry out a quantitative analysis of 
such a postulated system, it is essential to reduce the number of 
operational units to a minimum, otherwise the labor in carrying out 
mathematical calculations is excessive. However, the number of 
units selected is sufficient to demonstrate the basic phenomena, 
and provide a quantitative picture of the system. Calculations are 
carried out with certain reservations, but these are not arbitrary, 
and may well be in accord with certain experimental information 
available in regard to metabolic and synthetic processes. Chain 
reaction combined with synthesis can be expressed symbolically as 
the following, when a non-dissociable type of enzyme synthesis is 
used (Heinmets and Herschman, 1960): 


MaT SPP cteoes 


*NOTE. In the present analysis the end-product has no effect on en- 
zyme synthesis in preceding steps. Such may be a normal phenomenon. 
But it may represent a case where, as a result of the action of an ex- 
ternal agent, some intracellular enzymes have been destroyed, but en- 
zyme forming capacity is not affected; it may also be a case for a mutant 
cell where a block exists in enzyme synthesis and an intermediate sub- 
Strate accumulates. A recent study reveals (L. Gorini, 1960) ‘‘Antagonism 
between substrate and repressor in controlling the formation of a bio- 
synthetic enzyme,’’ Proc. Nat. Acad. Sci., 46, 682-690, that end product 


may affect enzyme formation in preceding steps. Such processes are cur- 
rently analyzed. 
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FIGURE 1. A scheme of sequential induction of enzyme synthesis 
coupled with linear chain reactions, 
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We will assume, for simplicity, that all the reactions occur only in 


one direction. Such an assumption, for metabolically and synthe- 
tically active cells can be considered a relatively accurate ap- 
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proximation (Krebs, 1947), The following notations will be used to 
indicate the concentrations of operational units: 
1. Repressed template (T)-a@ 4. Active template (TS;) ~ a; 
2. Pool (P) - 6 5. Enzyme (E,) - 0; 
3. Substrate (S;) -— s; 


Then concentration changes of the individual components can be 
expressed as follows: 


da 


dt = —k,s,;a-—k482a 


As say —k,g a,b —ksaob 


Gp 7 1810 -kaarb 


d 

TD = kasaa - ks aad 

db 

FE = 2a16 (1) 
db 

FT = ksaab 

LoL Eay 8,6, —kg82a—keg 82h 
ake on SST 4°2 

Bll ete 

ds, 

FR 


One immediately deduces the following conservation laws for the 
system: 


S; = constant 
b + by + bg = by = constant (2) 
@+@,+@2+ 6; + bg = ao = constant 


These give the connection between the initial and final state of the 
system. At the starting of the synthesis, the system consists of 
pools and initial substrate. The concentration of the latter is con- 
sidered to be constant during the synthesis—a condition which is 
conceivable in a metabolic system, i.e., an external substrate (in- 
ducer) concentration may be high enough and not significantly used 
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up in the synthetic process or the level of an internal substrate 
may be maintained constant by feed-back control mechanisms. 

For simplicity it is assumed that all rate constant values are the 
same. This assumption does not change the general validity of ex- 
pressions derived, but is rather a means to obtain equations which 
can be handled for numerical calculations. Furthermore it is as- 
sumed that the rate constants for the enzyme synthesis are much 
smaller than those of enzyme-substrate reactions. This is a na- 
tural condition and the subject has been discussed elsewhere 
(Heinmets and Herschman, 1960). The following relationships are 
indicated: 


ky ~k3~kg = ke 


Ill 
a 


(3) 
ko ks = yk 


where 
es oe be 


Concentrations of all individual units will be measured in relative 
terms of first substrate concentration. 


y=a/s, z= b/s, U2 = 82/81 
Y1 =@1/81 21 = bi/s1 ug = 83/81 (4) 
Yo = a2/81 29 = be/s1 
also put: 
= 40/81 B = 60/81 


and « =k s,t, where ¢ is time. 
In terms of these new variables, equation (1) becomes 


yo = —(1+ va)y 

2° = —y(yi + y2)2 

Yi=y-yyi2 

o =uU2y-—yy22 (5) 
Gee MEY aay 3:2 

2f=yyi2 

23 =yy22 


Ug = 21 — U2Y — U2Z2 


Us = U222 
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where differentiation with respect to z is indicated by prime. 


eral conditions are: 
2+21:+22=8 
YR¥ Tyg si 29 = Oo 


Initial conditions are: 


«= 0; 2= 8; y= a 
i = 2g = Yi = Y= tg = ty = 0 


The following interrelationships can be established: 


Yr1+y2= X-y — (21 + 22) 
= ny he 8) 
=a-Bt+2-y 


2°=-y2(a-B+2-y) 
n= ~y2(a-B-y)-y2? 


This equation has an immediate solution in terms of 


e()= [y+ Ba) ae 
0 


as 
Ber? 


x 
1+ By [ ae 
0 


2 = 


Where we have used condition 2 = f, 2 = 0. 
Also as: 

yityiy2=y 
So that, 


-{*y d x Y zdx 
vice ™ fyeb dz 
0 


but we have 


ao = —y8(0— 8 —y) ye 
or 


ye=y(y+B-a)-2’/2 


a ad = y — log 2/B 
0 


Gen- 


(6) 


(6°) 


(7) 


(8) 


(9) 


(10) 
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Consequently: 
i Y zdx ‘ B ve 
2 
“x Cie) 
=] + By | e” dz 
0 
and 


1 


v1 eee hry | oan) de (12) 
1+ By [ e% Pde o 9 
0 


Since 


if al aa) dem ["yde [9 de— ["(er# [” yas] dx 
0 fo) 0 0 0 0 
yon [vie ——_Br_[" of" yas) 
ety [ora 
0 


but 
[ yae= 9+ -B)2 
0 
Conseauently: 
Y= p+ (%— B)2 -——BY _ [eres (a - B)z) de 


1+ By [ ePtda 0 
; (13) 
Yi= P+ (A-B)a - yar” | eX? [p+ (a - B)aldz 
0 


Other components of the system can be derived using relationships 
from equations (5). 


ar+yr=y 

z (14) 

aityn [ yde=9+(%-B)@ 
0 


By substitution we find: 


2, = ya 7? [i eX? [p+ (a- B)aldz (15) 
0 
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22 =B-2- 41 


=p aftr yer? [a9 (9 «(a -praldel ce 


0 
Yg=X-B+2-y-y1 
=2-y¥i- 9 


Yo=—-(4- B)e- 9- y+ efaryer7? [e7Plo+ 0 - B)z] ie} 


0 


(17) 


To determine substrate concentrations we use again the previous 
relation (5) 
Ug = —(1+ y’/y) 


4 


Us = 21 —Uel(y + 22) = 21+ Y+ 22 + all) + 29) 


, 


B-2+yty t+ re 


I 


(18) 
va = Be - [ ede + [ ydesy-a+ | Zod (log y) 
0 0 0 
Bl Z me 
= 9(t)+ar-9+ —log>+y-A4+4+ 22d (log y) 
Nea: x 
vas a(e-teys gy tlog (1+ By [7% +s) + 
0 
aery [ e%P[p + (a - B)aldz (19) 
x , 
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3 222 a- ‘7 °72 
Ug = 21 — U2Y — U222Q 
(20) 


Ug + Us ~[ardes {ude + [yas 


af 2,dr+p+(A-B)r+y-a 
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An equation for y is then found by substituting 
Ug=-1- 7 into the equation (9). 


The resulting equation is a non-linear integro-differential equation, 
and cannot be solved in general. 

Before proceeding to find approximate solutions, we will explore 
the behavior of the individual components as function of 2 values. 
We note that in steady state (Asymptotic, i.e. as 7 — »), the de- 
rivatives of all components, except us’, become zero. 


y (1+ ue) =0 
2(y1 + y2) =0 
y=0; yiz=0; yo2=0 
21 = U222 
and 
per? 
1+ By [ ey? dx 
0 


se 


In conditions when y — 0 


~ — constant + (B - 4)z 


Let us put 7 =! ydex 
0 


It is assumed here, that y becomes zero rapidly enough for 7 to be 
finite. For approximate evaluations of various components of the 
system, we will distinguish three cases: 


1, = 6;22,-2.0 
Il. a> Bp; 2—-0 
Il, a<B; z2—~B-a 
Proof: 
O< psn 
lise? Sel 


x 
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As e”? is bounded when z — ~; 2 — 0 


Il. e%? —0 as e+ (P— -—-~) and [ eae 0 
) 


gia Ome aee es aco 


0 AP ss as 2— oo (P — +.) 


z2—y+B-a4 

z2— B - a, i.e., is bounded, as  — 
Since, 2+ 21 + 22 =B 

Y+Yit+y2=A-2;-22=A4-BHh+2 


Then for the assumed cases: 


I a=, z2—0 21+ 22 —-f8 
yD Yi+y2— 0 
Yi — 0 y2 — 0 
(Since neither can be negative) 
I. a8 z2—0 21+2,— 8 
tw Yitys— 4-6 
Ill. a<£B 2—B-4; 24+ 22 — & 


yO yity2—0; yi 0; y2 +0 


We also have the relation 
U2 — 21/29 


2 2 
hg ek) a 1___ case I, Il 
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Solutions for limiting cases can now be determined. 


Case I. 
yf e”? pda 
0 


x x 
1+ By [ oe dey [ ey? pdx - e%? 
0 0 


1 


As 9 is bounded, the integrals are unbounded. Using L’ Hospital’s 
rule, 


i ge laa Le Ah 2 
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(21 + 22) — 218 
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Case II. 
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ey? —+ 0 


All integrals are bounded. This case requires special treatment. 


ave Tene 


We note that it is also the case for which 


yi+y240, (e&—~). 
We note: 


2, + Yy1 — 7, in general and as 


ay tye =uay=-y-y 
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x 
Leos) ydxr+ua-y 
0 


22 + yg — 4-7, in general 


Accordingly the confusion in Zz; + 22 is the same asin yi + Yo- 

A general solution would yield information of how various opera- 
tional units of the system will be built up at various relative values 
of « and B. From the practical point of view this would mean that 
the relative size of template pool (nucleic acids) and synthetic 
pool (amino acids) would determine the transient-state formation- 
pattern of the enzymes and substrates. Such an analysis would be 
interesting, but perhaps at the present a little premature, since we 
do not have experimental data in hand, which would permit us to 
assign definite quantitative values to specific pools. Neverthe- 
less for the analysis of the system it is essential to derive definite 
solutions for the functions of the individual operational units. 
Since equations are easier to solve in the case where two pools are 
equal, and these solutions appear in a form which can be utilized 
for numerical calculations, we proceed with the case where 4 = B 
and y << 1, (as previously assumed). The latter assumes that en- 
zyme formation is very slow when compared with enzyme-substrate 
reaction, and this is in accord with experimental observations. 

We will now proceed to determine the approximate formation of 
the individual operational units of the system. 


We identify: Pp = ot =[ ydx (8’) 
0 


and, 


x? 
ene ee eee (9’) 


x 
1+ ay f eX? dx 
0 


ay [ pe”? da 
CP yn ia Sn (13’) 
eds ay f e”? da 
0 


Since the following relations exist: (on the basis of previous deri- 
vations) 
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Expressions for the individual units will be as following: 
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To proceed further, we note that, 
Py 
So that if, 
ya<<1l and yn<1 
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(These last are mutually consistent, as will be shown), we may put, 


eyP =1l+yQ 
Consequently, 

ie a(1+y) 

1+ aye 

ay | pdx 
212 pe eres 

: 1+ dye 

ay (x2 -9- f edz) 
0 
aa= 1+ay2 
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ay f pdxe—(1+ Aye) ue’ 
0 
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x 
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(23) 


(24) 


(25) 


(26) 


(27) 


(28) 


(29) 


(30) 


We are now in a position to prove a very important relation, namely 
that the concentration of intermediate substrate (S2) eventually be- 


comes equal to the concentration of the initial substrate (S;). 


Proof. From equation (P. 15, Case I) we have the relation (ug non- 


decreasing) 
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where [ P= eo yaa 
0 


By using equation (29) and noting if uy is slowly varying, we may 
set 


v - 1-€&; y= Qe ~ Ut) x (31) 
ene. Ge TT Ute) 5] 
? ‘ey a 
itp ee ges [1 — e748] 
IP SEO Seer aaa emersy 
if 
me ee 


Consequently approximate values for other components are: 


Bate Ay ar -P 
ae 1+€/\1l+aye 


=(f2) (sB =e 
1+€/\1+ayz2 


Hence by (5), we note 


and by conjecture in (31) 
ua € 
~€=1 and ug —~1=7/A-7 
So that 
7 = 4/2 [which justifies conjecture leading to equation (23)]. 


Corollary 
This means that after an initial rapid transient, the enzyme con- 


centrations become the same. It will be shown that it will occur 
in the order of the fast and not the slow reaction time. (It can be 
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shown that the previous corollary for enzymes is also true for ‘‘TS”’ 


complexes). 
On the basis of equation (31), we may then put (317) 


x ae 1 (24’) 
[ odes S39 
0 
21 [Zo a Ce — 0) 
2(1+ Ay x) 
Further relations are: 
~ P+ 7 ty (~- zy) 
y¥iz pa Bie ue EP (26) 


1+Qayz2 


The foregoing results make it possible to solve for ug to a higher 
approximation and offer some critical measure on validity of evalua- 
tion. Using (29) and (31’) 

ua = — (log y)” = 0 
So that from (28) we have 


a 
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Ug = Oo 
nd al mae 
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where the neglected terms are not significant. 


For an approximate analysis it is instructive to evaluate the con- 
centration of enzymes and substrates at three time zones: 


I. e < }- (time is in the range with the reaction 


T +s — TS) 


Il. 1 <a< nes (intermediate times) 
2 ay 


ie > = (time is long compared with the reactions: 
TS +P — E) 


The first region represents the rapid transients and here we have 
the following approximate relations: 
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1 2(1+ aya) Al 

a 

y (ae — #) 

m= Gl ~ Zilog (1+ aya)| - 9 


In the last zone we have: 


12 


U2 


r 


4 ae 
aha 
1+ ayz 
Y1 ~ ya ~ 2/2 
2 
a; “ite a a Ba 
2(1+ ayz) 2 
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bot Mag ee 


In order to illustrate the time dependent behavior of formation of 
substrates and synthesis of enzyme, an example will be used. De- 
rived equations for the region I, II, and III have been used to evalu- 
ate the individual components (for the arbitrary values « = 1 and 
y = 0.01), and some of these are presented in Fig. 2, as function of 
2 (equivalent to time). Components 2;, 22, wg and w3 are repre- 
sented as relative amounts to initial substrate concentration, which 
is taken as unity. Accumulation of end-substrate (us) is so rapid, 
that for this a different concentration scale has to be used. Forma- 
tion of substrate template complex is not graphed, since the process 
is completed very rapidly (z.e., for e = 3 y < 0.05) and does not fit 
into scale with other processes. A steady-state value is reached 
for intermediate substrate (w2) in a relatively short time. The in- 
crease of enzyme concentration does not change the concentration 
level of substrate and the increasing enzyme content produces an 
increasing flow of substrates via pathway. Formation of enzymes 
(2; and zg) proceeds almost simultaneously and will be completed 
when the pool is exhausted. After that, when the initial substrate 
is kept constant, a steady state flow of substrate proceeds via 
pathway. Rapid accumulation of end-substrate would produce ex- 
cessively high concentration of the latter. Such high concentra- 
tions can be inhibitory or be interfering with other metabolic reac- 
tions. However, the end-product could be used up continuously by 
some other reaction and so maintained at much lower concentration, 
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FIGURE 2. Sequential enzyme synthesis and substrate formation as a 
function of time. 


or it could be controlled by a feed-back mechanism (the latter will 
be analyzed in a subsequent publication). This sequential system 
contains only a few operational units, but in principal it has similar 
characteristics of a more extended system. Such a system could be 
a long chain of enzymes arranged in a sequential reaction series. 


Discussion and Summary 


It has been demonstrated that postulated schemes for sequential 
enzyme synthesis when coupled with linear chain reactions can be 
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used for mathematical calculations and these yield the information 
which describes enzyme synthesis and substrate transformation as 
a function of time. Analysis of such postulated systems may en- 
able us to understand and follow complex biological processes 
which cannot be understood solely on a descriptive basis. Quanti- 
tative transient-state characterization of enzyme synthesis coupled 
with enzyme-substrate reactions is our first step in the study of 
more complex cellular processes. These may be important not only 
for understanding normal growth but also for growth processes lead- 
ing to abnormalities. 

An analysis of linear chain reactions coupled with enzyme syn- 
thesis is interesting in itself—strictly from a phenomenological 
point of view—but in addition it may help to analyze certain ab- 
normal cellular processes. In the case of cellular injury, one can 
consider that various enzymes and substrates may have been in- 
activated or destroyed when an external agent interacts with meta- 
bolic system. One can ask, what will happen subsequently to 
these functional units when an injured cell is in a metabolically 
suitable environment where synthesis can take place? This prob- 
lem, which is very complex, contains many possibilities. For an 
example one can postulate that a build-up of enzyme systems may 
be initiated and may proceed as a sequential synthesis. Thus 
many functional processes which failed can be reestablished and 
the cell may acquire the ability to function again, but the system 
may be different from the original state. Which of the factors de- 
termine functional stability or alternative functional patterns can 
be analyzed at present only in postulated systems. Since we at- 
tempt basically to simulate actual systems, such studies may ‘‘open 
ways” to gain more insight and understanding in cellular phenomena. 


We would like to thank Dr. Selig Starr for carefully checking the 
mathematical solutions in the manuscript and for his suggestions. 
We are also grateful to Dr. Ernesto Trucco for reading the original 
draft of the manuscript and for his comments. 
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By introducing a plausible model for the initiation of axonal impulses 
the output is obtained as a function of the input incoming impulses. If 
the temporal aspects of the excitatory process resulting from the afferent 
impulses are sufficiently rapid one obtains the discontinuous or micro- 
scopic model of McCulloch-Pitts. If these are sufficiently slow a con- 
tinuous model, such as Rashevsky’s one or two factor theory, is a natural 
model. But the linear relation between the strength of excitation of one 
axon and excitatory factor of the next will not in general hold. However, 
under conditions which are not too restrictive the linear relation with 
threshold can be considered as satisfactory approximation over a fairly 
wide range of values. 


Recent work in neurophysiology has brought out a great deal of 
information on the mechanism of synaptic transmission (cf. Grund- 
fest, 1957; Tobias, 1959). It is the purpose of this note to con- 
sider how these findings affect the interpretation of mathematical 
models of neural net activity. The continuous two factor model 
was first introduced by N. Rashevsky (1948) as a formal model. 
While the discontinuous nature of impulses was considered, it was 
assumed that the excitatory or inhibitory states were continuous. 
Subsequently, in order to allow for certain observations the dis- 
continuous model (McCulloch and Pitts, 1943) was introduced in 
which output impulses followed, after a fixed delay, input im- 
pulses, provided that the latter were of sufficient number falling 
within a prescribed time interval. The requirements of quantized 
time and the uniformity of delay times could be relaxed under con- 
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ditions in which the impulses were uncorrelated in time (Landahl, 
McCulloch and Pitts, 19438). 

In order to reconcile the approaches a number of suggestions 
were made whereby the discontinuous model could be made to 
yield the results obtained on the basis of the continuous model. 
One of these (Landahl, 1945) pointed out certain similarities while 
others (Rashevsky 1945; 1948, Chap. 30) required that a macro- 
synapse be composed of many microsynapses together with postu- 
lated structural relationships. In view of recent observations 
such models may not be needed formally, although they may very 
well play a role in slower processes. 

In order to simplify matters we shall make the following assump- 
tions regarding the initiation of impulses in the axon. Let the 
potential of the axon near the cell body be V? when the axon is 
not excited. After excitation assume (1) that the potential re- 
covers from zero back to V? exponentially with a time constant T 
(Kamiya, 1958; Landahl, 1959). Let V; (¢) be the potential differ- 
ence across the cell body membrane at time ¢ and let R + Ry be 
the resistance as seen by the axon, where F, is the contribution 
due to the cell body membrane. Assume (2) that if the current 
density ¢ through the axon exceeds 7* the axon fires, its potential 
dropping suddenly to zero and then recovering toward V2. 

Regarding the properties of the cell body membrane and its 
afferent surfaces we shall consider three possibilities: I. It is 
electrically excitable, e.g. the membrane ‘‘fires’’ under certain 
conditions but the potential changes in the absence of excitation 
are unimportant. II. The cell body membrane polarizes due to in- 
coming impulses but does not excite. III. The cell membrane does 
not excite or polarize appreciably but impedance changes occur. 
We shall not consider combinations of these cases. 

In case I, it is likely that the time constants of the excitability 
of the cell body is rapid. If the processes are as fast as or faster 
than those of the incoming impulses there will be little or no tem- 
poral summation. At most, one output impulse will be produced by 
a Single incoming impulse or by several arriving nearly simultane- 
ously. If the delay times are sufficiently uniform and to the extent 
that one can quantize time, the model of W. S. McCulloch and W. 
Pitts (1943) can be used. If the delay times are not sufficiently 
uniform and the impulses are sufficiently uncorrelated in time the 


method of H. D. Landahl, W. S. McCulloch and W. Pitts (1943) can 
be used. 
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In case II it is more than likely that the important part of the 
polarization process takes substantially longer than the time con- 
stant T of the axon. If this were not so, however, one would use 
the same mathematical models as would be used for case I above 
if the process is fast enough. Hence we shall consider the case 
in which the polarization process is slow. In particular we re- 
Strict consideration to the case in which the potential V, changes 
little in the time T. Let the potential have some value V;. In the 
absence of incoming impulses we shall denote this value by V2. 
In terms of the model there will be at this time a current density 


t= (Va — Vo)/(R + Rp). (1) 


If this value of 7 is greater than 7* when V, = V2, a propagated im- 
pulse will be initiated. From the simplifying assumption intro- 
duced above we may write 


Ve =Va(i =e" */7), (2) 


The time elapsing between one impulse and the next is then given 
by solving for ¢ from (1) and (2) with z= 2*. If we denote this in- 
terval by ¢* we have 


t* = t log [V2/(V2 - i*R — i*Ry — Vp)). (3) 


In the above case we suppose that R, does not change even 
slowly so that we may introduce Rk’ =k + Ry. If we let 


1- Vi/Ve = €, (4) 


and suppose that V; > 0, then the output frequency will be given by 


1 a 
ee ee RV Oni OS hy. 
" lot log (6-4) +8 6) 


y= QOS fh, 


where A is a threshold value below which there is no response, 
and § is the time required for the fast process which we shall 
neglect hereafter. Note that / is the reciprocal of the ratio V2/R’ 
to i*, the maximum current to the threshold current, which ratio is 
essentially the safety factor for nerve. From (5) v increases from 
zero rapidly starting from € = A, then more slowly, but then again 
more rapidly near €=1. If one considers that the recovery when 
Vi is small is likely to be slowed down, this latter increase is not 
likely to occur. Furthermore, because of the finite recovery time 
which has been neglected, the frequency cannot continue to in- 
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crease. As a very rough approximation one could set (cf. Figure 1) 
= < 
ie ig | at OS (6) 
vy =A Toe Ay Ech, 


where A is a number of the order of magnitude of 3., Note, how- 
ever, that from relation (5) v rises rapidly near threshold. The re- 
lations (6) are just the expressions for the continuous model for 
neural elements introduced by Rashevsky to relate the frequency 
v with the excitation factor €. Hence in this case II we may at- 
tempt to identify € with the decrement in the polarization or the 
depolarization. 

There then remains the problem of interpreting the slower tem- 
poral aspects of € or of V}. We may think of the excitatory state 
as being the effect of a chemical mediator which is produced by the 
incoming impulses and which is eliminated at a rate proportional 
to the concentration at a given time. The rate constant for this 
process in this case would be the rate constant a in the continu- 
ous model (Rashevsky, 1945; 1948, Chap. 30) provided any possi- 
ble lag involved in moving from the production site to the site of 
action can be neglected in comparison to the time constant a}. 


vr 


E—> 


FIGURE 1, 
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On the other hand, if the action currents of the incoming impulses 
act directly on the post synaptic membrane producing changes 
which are restored with a rate constant a’, and if these changes 
result in potential (or impedance) changes with little or no delay 
(<< 1/a’), then one again would obtain the continuous model as a 
natural approximation. 

We consider next the case III in which the effect of incoming im- 
pulses is to produce transient changes in the impedance but no ap- 
preciable change in the polarization. If the time constants of the 
process are fast enough, then the results would be the same as in 
case I. But if the process is slow enough we may proceed as in 
case II. Let Vz be a constant which satisfies the relation 


i*R + i*R2>VO-Vy > G4R (7) 


where FP? is the value of Fy in the absence of incoming impulses. 
The above condition simply means that the potential difference is 
not so small that outgoing impulses cannot ever be produced or so 
large that they are produced even in absence of stimulation. It 
seems plausible to suppose that the conductivity of the membrane 
might change in proportion to the amount of a stimulating sub- 
stance. Hence let €’ be defined as the relative increase in con- 
ductivity: 
roeae 4 

Sal potas (8) 
If e°> fh’ = (*RP + *R -VP + V,)/(V2 - Ve - i*R) > 0, then there 
will be an output frequency v given by 


_- 71 dese Eda e VoeVe =k 
vie T Cea ail ae REA CRE aE (9) 
From the above expression it can be seen that the frequency rises 
slowly toward a finite upper limit. For small values of F, the fre- 
quency is relatively constant after the initial abrupt increase and 
hence the output frequency tends toward an all-or-more relation 
with the input. For F close to one (0< F <1), the frequency in- 
creases roughly linearly over a wide range after a small abrupt 
rise (Figure 2). 

If h’ is not too small, then a relation similar to (6) can be in- 
troduced as a rough approximation but with the coefficient A hav- 
ing, in this case, a value in the neighborhood of %. 

While it can be seen that under conditions which are not overly 
restrictive one can obtain a linear relation between the excitatory 
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factor € and the average output frequency, there remains a tend- 
ency for the output to show the initial abrupt rise near threshold. 
If one uses the linear approximation, then in a chain of neural ele- 
ments in the steady state, the threshold of the chain is the sum of 
the individual threshold weighted by combinations of amplification 
factors between elements while the slope of the input-output rela- 
tionship is proportional to the product of these factors. But if ex- 
pression (9) holds, the result is more complex. For example, for 
F =0.9, A’=1, the threshold of a chain of three identical elements 
in which each has the relation €(«) = ftv, where f= 83, is in- 
creased by a little over ten per cent and the output of the chain is 
rather similar to that of a single element over a wide range of 
values. But if the amplification factor, f, is decreased, the thresh- 
old of the chain decreases more rapidly than in the linear approxi- 
mation (Landahl, 19388). 

Another factor which has the tendency to make the linear rela- 
tion between € and v somewhat more plausible is the inherent 
variability in the input due to its discrete nature. The implied 
equation for the one-factor continuous model is 


Pie Ai(t) - ae (10) 
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where 7(¢) represents the effective stimulating current due to the 
incoming impulses. An average frequency is assumed which is 
high enough so that an average current, averaged over an interval 
of the order of 1/a, can be used. This average is 7 = /ov where [9 
is the charge moved per impulse, so that 7 = E in Rashevsky’s no- 
tation. If the number of afferent impulses arriving in the time 1/a 
is not large enough, there will be an appreciable fluctuation in the 
value of €, the relative variability being (va)-”, One effect of 
this variability is to modify the abrupt rise in the output frequency. 
If, for example, one replaces the variability by a rectangular dis- 
tribution covering the range 0.7 to 1.3 for the case of equation (6) 
with F = 0.9, kh’ =1, then the resulting output becomes almost a 
linear function of the excitatory factor € over a wide range of 
values (dotted curve in Figure 2). 
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Emphasis upon the importance of homeostatic feedback has drawn atten- 
tion away from the complexity of biological processes. A study of glucose 
metabolism indicates the importance of open-cycle as well as closed- 
cycle mechanisms. Besides the glucose-dependent mechanism of insulin 
secretion, many open-cycle processes involving the liver, adrenal glands 
and kidneys, play important roles in the variation of blood glucose. In 
addition, glucose utilization by the tissues is essentially open-cycle in 
nature. 


The process of homeostatic feedback strongly influences the be- 
havior of many physiological substances. A field of knowledge 
called cybernetics has sprung up to exploit the important role of 
feedback in many biological and social processes (Wiener, 1948). 
Stress produces, in this view, deviations from stable physiological 
levels and calls forth restorative forces. Changes in the values of 
the internal standards as the result of the stress may eventually 
damage the organism, in accordance with the concept of the ‘‘adap- 
tation syndrome’’ (Selye, 1950). Unfortunately, while homeostasis 
with feedback is a powerful concept theoretically and a useful 
concept practically in many instances, it frequently falls far short 
of the actual reality. The number of patterns besides feedback 
which affect the physiological variation of some substances is 
surprisingly large. Instead of dependence upon the closed-cycle 
system, this type of behavior is many-sided and could well be des- 
ignated by such a term as ‘‘mixed-cycle.’’ This paper analyzes 
the complexities of glucose metabolism in a mathematical sense 
as an illustration of the ‘‘mixed-cycle’’ concept. 

The blood sugar level ordinarily varies over a narrow range. In 
non-diabetics, it does not exceed 180 mgm/100 ml. At the other 
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extreme, the hypoglycemic level of 60 mgm/100 ml is infrequently 
reached. This stability exists despite the existence of irregular 
glucose loads and needs due to eating, exercise, and stress. Con- 
trol of the glucose level centers about the liver and pancreas. The 
release of insulin from the pancreas varies according to the level of 
blood glucose (White et al., 1959). A high glucose level leads to 
an increase in insulin secretion and vice versa. In addition, glucagon 
secretion helps to raise the blood level by increasing the amount of 
glucose released from the liver (Behrens and Bromer, 1958). The 
overall effect, mathematically, may approximate a restorative force 
proportional to the deviation of blood glucose from a standard value. 
Then, if g is the deviation of glucose from a certain standard value 


2 
G per 100 ml of blood, f —kg. In this situation, g would be an 


oscillatory function, with frequency /&/2z. In actual practice, oral 
intake of 1 gram glucose per kilogram of body weight in the glucose 
tolerance test results in an oscillatory variation of g. Normally, 
the peak is reached in approximately 60 minutes and the glucose 
returns to the initial value within 2 hours. Frequently, this is 
followed by an overshoot, with the glucose dropping below the 
fasting level for a variable period of time. Since the magnitude of 
the depression is much less than that of the maximum value, the 
curve resembles the motion of a spring with damping. Thus, the 
simple oscillatory pattern is inadequate for representation of glucose 
behavior even under carefully controlled conditions. The complexity 
of glucose metabolism is illustrated further by the fact that normally 
a second sugar load is more rapidly handled by the body than the 
first (Duncan, 1947). This increase in & is due, apparently, to in- 
creased secretion of insulin in response to the second load. Thus, 
k may be a function of g and ¢, instead of being a fixed constant. 

The role of the liver in glucose metabolism is quite complex and 
incompletely understood. Its removal leads to a progressive drop in 
blood sugar to lethal levels. Hepatic behavior has not been shown 
to be directly dependent, in distinctive fashion, upon blood glucose. 
Rather, the contributions of the liver to glucose metabolism depend 
upon such blood-borne substances as epinephrine, insulin, certain 
hormones of the adrenal cortex, and glucagon. The adrenal cortical 
hormones are particularly important, since the formation of glucose 
in the liver from amino-acid precursors—the process of glucone- 
genesis—depends upon such cortical steroids as hydrocortisone. 
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Indeed, the fasting blood sugar drops to hypoglycemic levels in 
adrenal cortical insufficiency. Now, the adrenal cortex 1ssto-a 
large extent, under the control of the hormones of the anterior pitu- 
itary gland. The output of hydrocortisone, in addition to its basic 
diurnal pattern of unknown origin, varies in accordance with the 
effect of stress upon the pituitary gland (Munson and Briggs, 1955). 
As long as it does not exceed limits which would be stressful to 
the organism, the level of blood glucose does not control the secre- 
tion of hydrocortisone. The behavior of the adrenal medullary hormone 
epinephrine is similar. The amount of this catechol amine normally 
present in the blood varies with glucose chiefly when low levels of 
glucose threaten the organism. The increased hepatic glycogenolysis 
which results from the increased presence of epinephrine serves to 
raise the blood sugar level in times of stress (Elmadjian e¢ al., 
1958). Thus, in a mathematical sense, since the blood level of 
glucose is under the influence of substances whose secretion it 
does not significantly control, a considerable portion of glucose 
metabolism is open-cycle in nature. This fact indicates an increase 
in the number of equations necessary for the mathematical deter- 
mination of blood glucose, to take into account the effects of in- 
dependent variables besides glucose and insulin. 

Another source of difficulty is the behavior of the kidneys, par- 
ticularly in diabetic individuals. No glucose ordinarily appears in 
the urine until the blood level exceeds 180 mgm per cent (White e¢ 
al., 1959). The urine of experimental animals given phloridzin 
contains glucose at much lower blood levels. Ordinarily, then, all 
of the glucose which enters the kidneys is reabsorbed and the effect 
of the kidneys upon overall glucose metabolism can be ignored. The 
appearance of glucose in diabetic urine in considerable amounts 
raises mathematical problems. The point at which glucose is first 
lost through the kidneys represents a mathematical discontinuity, 
since before and after it, the kidneys have a different mathematical 
effect. This moment must be determined for each situation specif- 
ically. The equations which apply in the absence of urinary loss 
provide a solution up to the appearance of glucose in the urine. Im- 
mediately after the beginning of glucose excretion, new equations 
which take renal loss into consideration must be set up, with initial 
conditions identical with the last values obtained from the previous 
solution. Each time glucose appears in, or disappears from, the 
urine, the same procedure must be followed. In this way, aadis- 
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continuity in organ behavior does not result in an irregularity in the 
solution for blood glucose. Although the kidneys do not participate 
in a closed-cycle system, they act in their own way to narrow the 
variation of blood glucose. 

The focus of attention thus far has been the blood level of glucose. 
The existence of glucose in the blood stream is, of course, only a 
means to an end. The end is the utilization of glucose in the body 
tissues as a source of energy. The body as a whole obtains its 
energy from a wide variety of substances. In fasting, fat is pref- 
erentially burned. Glucose is particularly important during intense 
activity when an oxygen deficit occurs, since the anerobic break- 
down of glucose to lactic acid releases energy in the absence of 
oxygen (White et al., 1959). The utilization of glucose is determined 
by the needs and preferences of the tissues of the body. These 
needs and preferences, however, are not, within a wide range, af- 
fected by the level of blood glucose. Rather, the quality and dura- 
tion of body activity is governed by infinitely varied combinations 
of voluntary desire and involuntary reaction which occur from mo- 
ment to moment. Only in clinical hypoglycemia and diabetic acidosis 
do the blood levels of glucose predominantly affect the behavior of 
the organism. Thus, with respect to glucose utilization, no closed- 
cycle system exists. 

Insulin facilitates the metabolism of glucose. It increases the 
entrance of glucose into the peripheral cell (Levine and Goldstein, 
1955). Evaluation of the overall role of insulin is complicated by 
the fact that nervous tissue does not require it and that its action 
in the liver is different from its peripheral effect. If g is the blood 
level of glucose per cc, it appears reasonable to assume, as a first- 
order approximation, that the rate of entrance of glucose into the 
cell will depend upon the value of g. It will, further, depend si- 
multaneously upon the amount of insulin fixed to the cell wall 
(Stadie, 1956), which may or may not be proportional at any given 
instant to the blood level of insulin. The rate of entrance of glucose 
into the cell may reasonably be represented by k(g — x)./, where & 
is a constant, / is the effective concentration of insulin within the 
cell wall and 2 is the concentration of glucose within the cell. 
While the release of insulin from the pancreas is proportional to 9s 
the effect of insulin upon a peripheral cell is not necessarily pro- 
portional to g, since such factors as the adhesiveness of insulin to 
the cell membrane must be considered. The utilization of glucose 
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by the cell, and therefore the value of 2 is also independent of g. 
Thus 2, g, and / are, in all probability, independent variables. The 


absence of a closed-cycle pattern in the utilization of glucose may 
lead to considerable mathematical complication. 
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SOUTHERN REGIONAL GRADUATE SUMMER SESSION IN 
STATISTICS AT VIRGINIA POLYTECHNIC INSTITUTE 


The 1961 session of the Southern Regional Graduate Summer Ses- 
sion in Statistics will be held at the Virginia Polytechnic Institute, 
Blacksburg, Virginia, from June 15 to July 22, 1961. Blacksburg is 
2100 feet high on the plateau in the Allegheny Mountains and offers 
a delightfully cool summer climate. 

The Virginia Polytechnic Institute, Oklahoma State University, 
North Carolina State College, and the University of Florida have 
agreed to operate a continuing program of graduate summer ses- 
sions in statistics to be held at each institution in rotation. The 
program was instituted at Virginia Polytechnic Institute in the sum- 
mer of 1954. 

The 1961 session, like previous sessions under this program, is 
intended to serve: 1) teachers of statistics and mathematics; 2) pro- 
fessional workers who want formal training in modern statistics; 
3) research and engineering personnel who want intensive instruc- 
tion in basic statistical concepts and modern statistical methodology; 
4) Public Health statisticians who wish to keep informed about ad- 
vanced specialized theory and methods; 5) prospective candidates 
for graduate degrees in statistics; and 6) graduate students in other 
fields who desire supporting work in statistics. 

The session will last six weeks and courses will carry five quarter 
hours of credit. Not more than two courses may be taken for credit 
at any one session. The summer work in statistics may be applied 
as residence credit at any of the cooperating institutions, as well 
as certain other universities, in partial fulfillment of the require- 
ments for:a graduate degree. The program may be entered at any 
session, and consecutive courses will follow in successive sum- 
mers so that it would be possible for a student to complete the 
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course work in statistics for a Master’s degree in three summers. 
Students must satisfy the remaining requirements for course work 
and thesis at the institution where they are to be admitted to can- 
didacy. The advanced courses may be accepted as part of the Ph.D. 
program of the participating institution. 

A limited number of fellowships will be available for applicants 
from certain specialized areas. Doctoral courtesy will be honored 
for those holding Ph.D. or M.D. degrees. 

The courses to be offered in statistics in 1961 at the Virginia 
Polytechnic Institute are as follows: Statistical Methods, Sampling 
Theory; Applied Statistics for Engineers and Physical Sciences; 
Theory I, Probability; Theory II, Statistical Inference; Theory III, 
Theory of Linear Hypotheses; Non-parametric Methods; and Multi- 
variate Methods. 

A number of courses in advanced mathematics will be available 
during the Summer Session. For a complete listing please consult 
the University timetable. A series of Colloquia involving recent 
developments in statistical theory and methods will be conducted 
during the special Summer Session. 

Requests for application blanks for the summer session and for 
fellowships should be addressed to: 


Dr. Boyd Harshbarger, Head 
Department of Statistics 
Virginia Polytechnic Institute 
Blacksburg, Virginia 
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SYMPOSIUM ON MATHEMATICAL THEORIES OF 
BIOLOGICAL PHENOMENA 


A Symposium on Mathematical Theories of Biological Phenomena 
will be held at the Barbizon-Plaza Hotel in New York City on May 
8, 9, and 10, 1961 under the joint sponsorship of the Committee on 
Mathematical Biology of the University of Chicago and of the New 
York Academy of Sciences, under the general chairmanship of N. 
Rashevsky. 

One hour technical lectures will be given on different topics of 
mathematical biology with ample time for discussions. 


List of Participants and Topics: 


A. F. Bartholomay, Harvard Medical School—Mathematical Theory 
of Reaction Rates. 

G. Segre, University of Torino, Italy—Applications of Transfer 
Functions in Pharmacology. 

J. G. Defares, University of Leiden, Holland—Theory of Lung 


Functions. 
S. Roston, University of Louisville Medical School—Blood Pres- 


sure and the Cardiovascular System. 
E. H. Kerner, University of Buffalo—Gibbs Ensembles and Bio- 


logical Ensembles. 
G. Sacher and E. Trucco, Argonne National Laboratory—Sto- 


chastic Theory of Mortality. 
D. A. Rigas, University of Oregon Medical School—Theoretical 
Considerations of Cellular Proliferation as Studied by Cell Label- 


ing Processes. 
G. Karreman, Eastern Pennsylvania Psychiatric Institute—Elec- 


tronic Aspects of Quantum Biology. 
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H. D. Landahl, The University of Chicago—Mathematical Theory 


of the Central Nervous System. 
W. Rall, National Institutes of Health—Theory of Physiological 
Properties of Dendrites. 
J. H. Woodger, University of London—Axiomatization of Biology. 
N. Rashevsky, The University of Chicago—Mathematical Founda- 


tions of General Biology. 


Session Chairmen: 


J. Newton Kugelmass, Consultant, Departments of Health and Hos- 
pitals, New York City. Herman Branson, Howard University. J. H. 
Woodger, University of London. John Z. Hearon, National Institutes 
of Health. H. Quastler, Brookhaven Nationa! Laboratory. H. L. 
Lucas, State College Section, Raleigh, N. C. 


The symposium is made possible by a grant from the National 
Science Foundation to the University of Chicago. 


THE MEETINGS ARE OPEN TO ALL INTERESTED. 


